UFDC Home  myUFDC Home  Help 



Full Text  
PAGE 1 1 A CURVATURE BASED METHOD FOR SYSTEMATIC ERROR S ADJUSTMENT IN AIRBORNE LASER SCANNING DATA By PRAVESH KUMARI A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT OF THE REQUIRE MENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY UNIVERSITY OF FLORIDA 201 1 PAGE 2 2 201 1 P r a vesh Kuma ri PAGE 3 3 To My Father PAGE 4 4 ACKNOWLEDGMENTS I would like to acknowledge the support of my mentors, peers, and staff in the Civil and Coastal E ngineering at the University of Florida. In particular, I want to thank Dr. Clint Slatton, Dr. Ramesh Shrestha, and Dr. Bill Carter for their candor and guidance througho ut my career in Geosensing System Engineering ; Dr. Jasmeet Judge Dr. David Bloomquist for guid ing me through the process, Heezin Lee, Bot, Bidhya Yadav Mike Starek Abhinav Singhania and Juan Fernandez for their support and inspiration; Sailesh and Mansi for their guidance and encouragement; and Nancy Been for registering me in cours es and helping with all that paper work I truly would like to thank you all. I would also like to thank John Caceres for the opportunity to work together on research and writing projects at the university. I would like to especially thank Dr. Carter for his time, effor t, and graciousness; Thelma and Laura for their support and to the online learning networks and communities in which I have participated and have yet to participate in. Finally, I would like to acknowledge my parents and my uncle for inspiring and trusting me all the way Priyanka Parth a and Prasant I thank you for making all of this possible. PAGE 5 5 TABLE OF CONTENTS P age ACKNOWLEDGMENTS ................................ ................................ ................................ .. 4 LIST OF FIGURES ................................ ................................ ................................ .......... 8 ABSTRACT ................................ ................................ ................................ ................... 12 CHAPTER 1 INTRODUCTION ................................ ................................ ................................ .... 14 Airborne Laser Swath Mapp ing ................................ ................................ ............... 14 LiDAR Data Accuracy and Sources of Error ................................ ........................... 16 Problem Statement ................................ ................................ ................................ 19 Overview ................................ ................................ ................................ ................. 20 2 BACKGROUND AND LITERATURE REVIEW ................................ ....................... 23 Definition of Accuracy for LiDAR Data ................................ ................................ .... 30 Statistical Accuracy Assessment ................................ ................................ ............ 31 3 LINEAR AND PLANAR FEATURES BASED STRIP MATCHING .......................... 34 Geore ferencing LiDAR Data ................................ ................................ ................... 34 Methodology ................................ ................................ ................................ ........... 37 Formulation of Error Estimation Model ................................ ............................. 38 Least Square Adjustment of LiDAR Strips ................................ ........................ 43 Sampling and Classification ................................ ................................ ............. 44 PDF Estimation with Parzen Window ................................ ............................... 45 Sample Collection ................................ ................................ ............................ 47 Boundary Extraction and Tracing ................................ ................................ ............ 50 Regulariz ation and Least Square Adjustment ................................ ......................... 52 Evaluation and Results ................................ ................................ ........................... 53 4 DATA SAMPLING AND FEATURE EXTRACTION FOR SURFACE MATCHING .. 58 Feature Selection ................................ ................................ ................................ .... 60 Polynomial Fitting and Curvature Estimation ................................ .......................... 61 Surface Matching ................................ ................................ ................................ .... 63 Results ................................ ................................ ................................ .................... 64 Conclusion ................................ ................................ ................................ .............. 69 5 CALIBRATION OF LIDAR DATA USING CURVATURE BASED SURFACE MATCHING ................................ ................................ ................................ ............. 71 PAGE 6 6 Proposed Method ................................ ................................ ................................ ... 74 LiDAR Equation ................................ ................................ ................................ 75 Surface Matching ................................ ................................ ............................. 76 Curvature Based Sample Selection ................................ ................................ .. 83 Normal Estimation ................................ ................................ ............................ 84 Outlier Removal ................................ ................................ ................................ 84 Computational Aspects of Weighted Least Square Matching ........................... 85 Convergence ................................ ................................ ................................ .... 88 Results and Discussion ................................ ................................ ........................... 88 Conclusion ................................ ................................ ................................ .............. 99 6 DISSERTATION CONCL USION AND FUTURE WORK ................................ ...... 101 APPENDIX SOURCES OF ERRORS AND THEIR EFFECT ON LIDAR DATA ...... 103 Beam Divergence ................................ ................................ ................................ 103 Target Geometry and Reflectance ................................ ................................ ........ 104 Platform Position and Attitude ................................ ................................ ............... 105 GPS Posit ioning ................................ ................................ ............................. 105 Exterior Orientation from Inertial Measurement Unit ................................ ...... 106 Bore Sight Angles ................................ ................................ ................................ 107 Scan Angle ................................ ................................ ................................ ........... 107 Laser Range ................................ ................................ ................................ ......... 107 Atmosphere ................................ ................................ ................................ .......... 108 Interpolation and Coordinate System Conversion ................................ ................. 110 Effects of Major Systematic Errors ................................ ................................ ........ 111 Attitude Errors ................................ ................................ ................................ 111 Roll error ................................ ................................ ................................ .. 111 Pitch error ................................ ................................ ................................ 112 Heading error ................................ ................................ ........................... 113 Scan Angle Errors ................................ ................................ .......................... 115 Range Error ................................ ................................ ................................ .... 116 LIST OF REFERENCES ................................ ................................ ............................. 119 BIOGRAPHICAL SKETCH ................................ ................................ .......................... 123 PAGE 7 7 LIST OF TABLES Table P age 1 1 Magnitude of horizontal and vertical error in reconstructed L iDAR surface due to various systematic errors; flying height 1000m (above ground) and maximum scan angle 20 degrees. ................................ ................................ ...... 18 2 1 ) ............. 33 3 1 Systematic error estimates for Hogtown dataset ................................ ................ 54 3 2 Comparison of correlation analyses of flat planar features only a nd bounded surfaces ................................ ................................ ................................ .............. 56 4 1 RMSE (in meters) for different terrain types with three different sampling methods ................................ ................................ ................................ .............. 67 5 1 Calibratio n parameters for lava surface data ................................ ...................... 93 5 2 Correlation matrix of the estimated parameters (left) and parameter values with their respective standard deviations (right) (meteor crater data set) ........... 95 5 3 Stability analysis of proposed method by comparing intentionally added (black) and recovered systematic errors [curvature based method in (red); commercial software in (blue)] (Study Area 2) ................................ .................... 97 PAGE 8 8 LIST OF FIGURES Figure P age 1 1 Airborne laser swath mapping (Carter et al., 2007) ................................ ............ 14 1 2 Distribution of systematic errors across a flight strip: scanner scale; attitude errors such as errors in roll, pitch and heading; range bias; cumulative effect of all the errors ................................ ................................ ................................ .... 17 1 3 (a) Artifacts in DEM due to systematic discrepancies (b) no artifacts in DEM after systematic errors adjustment ................................ ................................ ...... 19 2 1 Accuracy and precision ................................ ................................ ...................... 23 2 2 Distribution of elevation and positional discrepancies across the overlapping region between two adjacent strips on flat surface (top) and sloped surface (bottom) ................................ ................................ ................................ .............. 27 2 3 Comparison of error distribution over an overlap region (60% overlap) between two adjacent strips flown in (top) same and (b ottom) opposite directions: (left) across flight profiles of the overlap region; (center) top views of the overlap area; black scan line is ground truth; (right) relative errors between the overlapping strips: along flight position error (X error in black) across flight position error (Y error in red) and elevation error (Z error in blue) .. 29 2 4 Horizontal accuracy assessment s ................................ ................................ ...... 32 3 1 Reference frames for ALSM system (adapted from Schenk, 2001) .................... 34 3 2 (a) WGS reference ellipsoid in relation with local level reference frame (ENU), (b) vector relationship among different sensor frames in WGS reference system ................................ ................................ ................................ 35 3 3 Tie features ................................ ................................ ................................ ......... 38 3 4 Incidence angle of the laser beam and surface normal of the laser footprint ...... 42 3 5 PDFs o f mean of first stop elevation and surface roughness of non ground points ................................ ................................ ................................ .................. 46 3 6 Study area for accuracy assessment: classified LiDAR points based on surface roughness (brown: higher values of roughness like forest; green: lower values of roughness like planar rooftops) ................................ .................. 46 3 7 Extracted planar features: (a) extracted buildings (b) gable roof (c) cross section of gable roof (d) complex gable roof (colors indicate different flight lines) ................................ ................................ ................................ ................... 48 PAGE 9 9 3 8 Extracted features: (a) a road classified based on elevation (b) small road patch shown in two overlapping flight lines. ................................ ........................ 48 3 9 Extracted building with planar surfaces ................................ .............................. 50 3 10 Extracted building with planar surfaces ................................ .............................. 51 3 11 Modified convex hull (Sampath and Shan, 2007) ................................ ............... 51 3 12 Traced boundaries of a building roof ................................ ................................ .. 52 3 13 LiDAR data classified based on elevations (University of Flor ida campus) ........ 55 3 14 Extracted boundaries (University of Florida campus) ................................ ......... 56 4 1 (a) Six degrees of freedom in matching overlapping strips; (b) two ring neighborhood ................................ ................................ ................................ ...... 60 4 2 Performance of sampling algorithm on synthetic LiDAR data (fractal); valley points (blue) and ridge points (red). ................................ ................................ .... 65 4 3 Results of curvature based sampling on various terrain types: (a) sand dunes, (b) flat terrain with few features and (c) mountainous terrain .................. 67 4 4 Convergence analysis for three different surface types ................................ ...... 68 5 1 LiDAR data strip illustrating non uniform point density along the flight direction; an example of data collected with high flight dynamics. ...................... 74 5 2 LiDAR data strip with uniform point density along the flight direction implying low flight dynamics. ................................ ................................ ............................ 74 5 3 Sampling based on curvature change (a) DEM of a hilly surface; (b) all classified ridge points, and (c) down sampled ridge points (red) ........................ 83 5 4 Work flow of the proposed surface matching and systematic error adjustment algorithm. ................................ ................................ ................................ ............ 87 5 5 Three different types of terrain: (left) mountainous lava surface; (center) sand dunes; and (right) flat surface with few features. ................................ ................ 90 5 6 Study area 1: (a) intensity image (b) flight swaths ................................ .............. 91 5 7 Study area 2; (a) DEM of a strip; (b) LiDAR points were selected using curvature based sampling from a small area (red box on DEM); (c) overlapping strips with ground control (line across the strips) ........................... 92 PAGE 10 10 5 8 Elevation discrepancies over first overlap area. (left) estimated with commercial software (average elevation error = 0.060m) ; (right) with proposed calibration algorithm (average elevation error = 0.0012m) ................. 94 5 9 Elevation discrepancies over second overlap area. (left) estimated with commercial software (average elevation error = 0.11m) ; (right) with proposed calibration algorithm (average elevation error = 0.0026m) ................................ 94 5 10 Calibration results for three overlapping flight lines: (a), (b) and (c) are the histo grams of elevation differences between three strips and the ground truth before adjustment. (d), (e) and (f) are the histograms of elevation differences of same three strips compared to ground truth after systematic error corrections were applied to them. ................................ ................................ ....... 95 5 11 Convergence analysis of proposed method over three different types of terrains ................................ ................................ ................................ ............... 99 A 1 Transmitted and received laser pulse (CFD) ................................ .................... 104 A 2 Laser foot print on ground ................................ ................................ ................ 105 A 3 GPS positioning (a) absolute (b) differential (c) IMU gyros and accelerometers mounted orthogonally (Colomina, 2000) ................................ 106 A 4 Profiles of temperature and pressure with altitude ................................ ............ 110 A 5 resulting vertical (dz) and horizontal (dy) error ................................ ................................ ................................ .................. 111 A 6 splacement (dx) error ................................ ................................ ................................ ........... 112 A 7 Heading or flight line showing horizontal (dx and dy) errors ................................ ................ 114 A 8 Scanner mirror error ................................ ................................ ......................... 115 A 9 Scanner sc ale ( sc) error ................................ ................................ ................. 115 A 10 Range bias ( ) ................................ ................................ ................................ 116 A 11 Across flight line profile of two overlapping strips over a gable roof ................. 116 A 12 E ffects of systematic attitude biases (roll, pitch and heading) on the nominal surface (blue) are visible in the modified surface (red) ................................ ..... 117 PAGE 11 11 A 13 Effects of systematic scan angle offset and scale bias on the nominal surface (blue) are visible in the modified surface (red) ................................ .................. 118 A 14 Effects of systematic range bias on the nominal surface (blue) are visible in the modified surface (red) ................................ ................................ ................. 118 PAGE 12 12 Abstract of Dissertation Presented to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for the Degree of Doctor of Philosophy A CURVATURE BASED METHOD FOR SYSTEMATIC ERRORS ADJUSTMENT IN AIRBORNE LASER SCANNING DATA By Pravesh Kumari MAY 201 1 Chair: Ramesh L. Shrestha Major: Civi l Engineering To produce research quality airborne laser scanning (ALS) data, it is necessary to remove systemic errors using a comprehensive error adjustment model based on the knowledge of the specific ALS system, including details of the data flow within all the components and the possible systematic and random errors associated with each component This dissertation explains the nature of various systematic errors and their individual and cumulative effects on the final reconstructed surfaces obtained from the ALS data. Then a mathematic al model is formulated to incorporate all the observations from vital components such as the global positioning system (GPS), inertial navigation system (INS) and the laser scanning system. With this model as the basis, two adjustment methods are developed to regulate the systematic errors in the data. Both of the proposed adjustment methods take into account the terrain characteristics and availability of the suitable features necessary to achieve a well conditioned adjustment solution. The first method utilizes linear and planar geometric features, abundant in the urban scenes, to match overlapping ALS strips and compute PAGE 13 13 adjustment s to param eters to remove systematic error s. N atural terrains were typically found to lack the requisite linear and planar geometric features required for the first met hod to work well, and a second method based on natural surface characteristics was developed. Surface matching is a well researched topic in Compute r Vision as well as terrestrial LiDAR where the size of range images is much smaller compared to airborne L iDAR data. Iteratively closest point (ICP) and its variants have been successfully used to align and register multiple overlapping views of an object. However, this technique is still new for airborne laser scanning due to various implementation issues In this dissertation an effort has been made to address these issues and full automation of the adjustment algorithm has been i nvestigated. In addition a technique has been presented to adjust overlapping strips using geometrical attributes such as curvatu re changes in a given terrain The magnitude and direction of these curvature changes in the data points indicate their ability to constrain the relative movement between the overlapping laser strips. The points from overlapping strips are matched through modified point to plane based iteratively closest point (ICP) method The suitability of the curvature based surface matching method has been tested for various terrains such as hilly, flat and rolling terrains. Finally, the dissertation concludes that a f eature based method can be successfully used in estimating the systematic biases present in the ALS data collected over natural or manmade terrain. PAGE 14 14 CHAPTER 1 INTRODUCTION Airborne L aser Swath M apping Airborne Laser Swath Mapping (ALSM) also known as Light Detection and Ranging ( LiDAR ) is a fairly new technique in the field of active remote sensing for acquiring more comprehensive and precise topographic information than traditional methods. The method relies on measuring the distance from a platform such as an airplane or helicopter trip travel time o f a brie f pulse of laser light ( Figure 1 1 ) The travel time is measured from the time the laser pulse is fired to the time laser light is reflected back from the surface. The reflected laser light is received using a small telescope that focuses any collected laser light onto a detector. The travel time is converted to distance from the plane to the surface based on the speed of light. Typically, a laser tr ans mitter is used that produces a near infrared laser pu lse that is invisible to humans (Airborne 1 Corporation, 2001). Figure 1 1 Airborne laser swath mapping (Carter et al., 2007) PAGE 15 15 LiDAR technology integrated with the Global Positioning System ( GPS ) and Inertial Measuring System or Un it (IMU) are mounted on an aircraft and flown over the study area. Laser transmitters are used th at fire tens to hundreds of thousands of pulses per second and record both the round trip travel time for each pulse and the a ngle at which it is reflected. By scanning the laser pulses across the terrain using a n oscillating mirror, a dense set of distances to the surface is measured ( Figure 1 1 ) along a narrow swath The LiDAR observations are converted to map coordinates and e levations for ea ch laser pulse by combining vector range data i.e. the laser ranges along with scan angles, with information on the position and orientation of the airplane at the time the laser pulse was fired. GPS provides positional r eference to a cert ain coordinate system such as WGS84 and IMU measure s the attitude angles i.e. roll, pitch, and yaw of the aircraft. The major advantage of ALS is that it directly gives the elevation of the scanned object points. D etailed height profiles of the vegetation can be created and because s ome of the laser light makes it through the openings in the foliage and reach es the ground. I t is also possible to get the coordinates of ground points under vegetation. The number of laser pulses reaching the ground largely de pends on the foliage density. Once the ground returns are identified, they are used to produce what is known as a digital elevation model (DEM) that describes the ground topography using a regularly spaced grid of elevation values. The data can also be use d to determine the height and density of the overlying vegetation, and to characterize the location, shape, and height of buildings and other manmade structures. Digital aerial photographs can also be taken while LiDAR is flown, provid ing an additional lay er of data PAGE 16 16 T he LiDAR data used to test the proposed algorithms was collected with an Optech Inc. Gemini unit The Gemini system owned and operated by National Center for Airborne Laser Mapping (NCALM) is capab le of collecting four returns per laser s hot, at two different settings of the divergence of the laser beam i.e. 0.25 mrad and 0.8 mrad The pulse repetition frequency ( P RF) can be set at 33 75, 100, 125, 140 or 167 KHz. LiDAR Data Accuracy and Sources of Error The accuracy and reliability of LiDAR s ystem s has become an active area of research, as many applications require decimeter level accuracy in topographic data sets. The assessment of systematic errors present in LiDAR data is a complex process because there are several sources of error t hat contribute to the gross error. A part from the target reflectivity properties and laser beam incidence angle, the main limiting factors are the accuracy of the platform position and orientation derived from the carrier phase differential GPS/INS data, s canner pointing errors and uncompensated effects in system calibration. The error budget of the mapping system can b e expressed as the contribution of the errors from core sub systems; laser rangefinder, IMU or Inertial Navigation System (INS) attitude sol ut ion, GPS position solution, latency of the receiver electronics, sensor mounting bia ses, scanning angle errors etc. ( Airborne 1 Corporation, 2001 ) Certain types of errors such as GPS offset and roll angle error simply translate or rotate the laser po int cloud, whereas errors such as a scan angle error, range bias and pitch angle bias tend to depend on the variations in other parameters ( Filin, 2003 ) Therefore, the overall accuracy and precision of the data obtained from the system depends on the ass embly of its different components. PAGE 17 17 For the above mentioned reasons, ALS observations are usually corrupted with systematic and random errors to some extent (Shrestha et al., 1999). Schenk (2001) and Triglav Cekada et al. (2009) have explained the effect of the systematic errors on positional and elevation accuracy of LiDAR ground points in detail. Figure 1 2 illustrate s the effect of five types of systematic errors in position and elevation of a reconstructed LiDAR surface across the flight direction. The se errors are scanner scale bias, attitude errors and range bias. It is clear that the result of scanner scale or range bias in the system is a non linear warped surface across the data strip. In addition, the magnitude of positional discrepancies is great er than the errors in the elevation. Figure 1 2 Distribution of systematic errors across a flight strip: scanner scale; attitude errors such as errors in roll, pitch and heading; range bias; cumulative effect of all the errors The horizontal as well as the vertical accuracy of L i DAR data depends on flying height. A vertical accuracy on the order of 15 cm root mean squared error (RMSE) has PAGE 18 18 been reported by most of the LiDAR vendors (ASPRS LiDAR Guidelines 2004) If flight conditions are optim um vertical accuracies of 7 8 cm RMSE have been reported (Lee et al., 2007). The a ctual accuracy of the LiDAR point cloud s and the interpolated grids or digital elevation models (DEMs) depends on many factors Therefore different studies report the vertical accuracy of LiDAR data with different results. For most of the studies the LiDAR data a re collected under leaf off conditions (Shrestha et al. 2001). The accuracy of the data collected in leaf on conditions, especially for the vegetation studies, is lower compared to data collected in leaf off season. The n umber of laser pulses that penetrate the vegetation during the leaf on season is reduced. As a result less number of pulses are reflected back from the ground, which significantly reduce the accuracy of the LiDAR data and the derived DEMs. Table 1 1 Magnitude of horizontal and vertical error in reconstructed LiDAR surface due to various systematic errors; flying height 1000m (above ground) and maximum scan angle 20 degrees. Magnitude of Systematic Error Maximum H orizontal Error (m) Maximum Vertical Error (m) Roll Angle (0.01 degree) 0.17 0.06 Pitch Angle (0.01 degree) 0.17 0. 0 6 Heading Angle (0.01 degree) 0.17 Scanner Scale (0.00 0 1) 0. 1 0 0.1 0 Laser Range (0.01 m) 0.01 with no minal aircraft attitude ( no flight dynamics ) The variations in the accuracies reported in different studies can be attributed, in part, to the differences between laser systems employed, flight characteristics, the terrain surveyed, how well the LiDAR was able to penetrate vegetat ion, and differences in the processing of the data itself such as differences in the algorithms used to filter PAGE 19 19 out returns from vegetation In short the vertical a ccuracy ranged from 3 to 100 cm, with the majority of the studies reporting from 7 to 22 cm (ASPRS LiDAR Guidelines 2004 ; Latypov 200 2 ) Proble m Statement As discussed earlier, the LiDAR system observations are often corrupted by systematic and random errors. Most of these errors are kept under control by implementing various strategies at t he flight mission level. To minimize the atmospheric effects, the aircraft are flown at relatively low elevations usually below 1000 m above ground level. To keep the GPS and INS errors in control, the missions are flown in short strips and during periods of low flight dynamics. However, it is not always possible to avoid such error causing factors and when present they shift, rotate or scale the actual position s of the laser point s on the ground. (a) (b) Figure 1 3 (a) Artifacts in DEM due to sy stematic discrepancies (b) no artifacts in DEM after systematic errors adjustment PAGE 20 20 The effect s of these errors are visible in the reconstructed surfaces or DEMs derived from the LiDAR data. Since, the LiDAR data is collected in long strips with as much as 50 60% overlap, there are di screpancies in positions and elevations of the conjugate feature points belonging to the same target in the overlapping strips. The physical artifacts such as lines or patterns are visible in the reconstructed surfaces in overl ap regions (Figure 1 2). The main objective of this dissertation is to fix these discrepancies at the system level and process ing the data free of systematic error s This involves understanding the effect of various systematic errors on the reconstructed s urface, proposing an adjustment model and testing its performance on different types of surfaces. In subsequent chapters, the adjustment procedures are explained in detail. Overview The objective of this dissertation is to identify major sources of errors in the airborne laser scanning (ALSM) data and propose a comprehensive adjustment technique to compute the systematic error estimates or calibration parameters during post processing. Chapter 1 introduces the concept of an airborne laser scanning system a nd its various components. It explains the motivation behin d this study and the need of an error adjustment algorithm that can comprehensively estimate the systematic errors in the data. Chapter 2 elaborates on the background of the error adjustment method s and narrates the previous work that has been done in the area. There has always been a need to establish a standard system to define the vertical and planimetric accuracy of the LiDAR data. However, even with point densities of several points per square meter achieved by the newer airborne LiDAR systems, the chances of two LiDAR points hitting the exact same spot on the ground are minimal which provides little redundancy in the data points. It is demonstrated in Chapte r 3 that by using linear PAGE 21 21 and planar surfaces rather than only planar surfaces in the urban scenes, it is possible to obtain good estimates of the systematic biases in the system measurements Th e basic equation for the geo referencing of the laser point to the suitable local coordinate fra me is explained in this chapter and a different approach from the existing and widely used calibration method s is proposed. It is concluded that the use of breaklines, along with the planar surfaces helps with the reliability of estimated calibration param eters. In addition, Chapter 3 explains the steps involved in implement ation of the proposed technique which has been successfully tested on LIDAR data collected over the University of Florida Campus This chapter also emphasizes the importance of robust classification and feature extraction methods. The proposed method, in this chapter, is specific to the urban topographic data and may not be useful for a more natural environment. Chapter 4 lays the foundation for another adjustment method that can be u sed for natural terrain. made features as the previously proposed method does. The importance of curvature based features present in the natural terrains has been emphasized for adjusting the systematic biases in the LiDAR system A performance based comparison of the curvature based sampling method has been made for the two most popular sampling techniques i.e. uniform and random sampling. In Chapter 5, a systematic error adjustment model based on curvature based sampling is propos ed. The successful application of the proposed algorithm is demonstrated by adjusting the systematic errors in different types of terrains including flat terrain, rolling surface and mountainous terrain. Finally, the results of this stud y are summarized in Chapter 6 and the foundation is laid out for future work. Appendix A PAGE 22 22 explains the sources of different types of errors and illustrates their effect s on final ground coordinates derived from LiDAR data. PAGE 23 23 CHAPTER 2 BACKGROUND AND L ITERATURE REVIEW As th e number of applications using LiDAR data has in creas ed, the accuracy of LiDAR has become an increasingly important issue for researchers ( Shrestha et al 1999 ; Airborne 1, 2001 ) The theoretical accuracy of a LiDAR system can be computed by incorporating the precision of each individual compone n t of the system. However, for a number of reasons, the accuracy obtained in the field is usually wors e than the specified theoretical accuracy (Latypov, 2002) Since t here has been lack of a common accuracy standar d for LiDAR data different vendors and users measure and report the overall accuracy of the collected data differently. Latypov (2004) mentioned that several times the confidence level (c.l.) associated with the accuracy measurement is not specified. The National (USA) Standard for Spatial Data Accuracy (NSSDA) defines accuracy a s being measurable at 95% c.l. and assumes errors have a normal distribution. In real life, the assumption of normal distribution is often rather poorly realized Figure 2 1 Acc uracy and p recision PAGE 24 24 Sometimes the accuracy associated with 95% c.l. is assumed to be double of the 68% c.l. accuracy to meet the NSSDA accuracy specifications. However, this is not a correct practice given the distribution of errors is normal or bell shape d. Accuracy or reference ( Figure 2 1 ) The difference between tests, by an order of magnitude than that of the one being tested. The root mean square erred to as absolute accuracy. Errors can be characterize d in three main categories, including blunders, systematic errors, and random errors. Blunders are significantly larger in magnitude compared to t he other two types. They can be detected and elimina ted with use of redundant observations. Systematic errors are caused by imperfect instruments, or deficiencies in the mathematical model used to compute the desired parameters from the observations. With appropriate surveying methods, they can either be el iminated during the observation process, or determined and corrected. This requires redundant observations, and to some extent, independent control information Random errors are always present and can never be eliminated. They can, however, be minimized, again by redundant observations. provides a measure of the theoretical accuracy, also referred to as precision. For the case where the estimated errors are purely random, i.e. neither blunders nor systematic PAGE 25 25 errors exist ; the the oretical accuracy and the absolu te accuracy are equivalent ( Friess, 2006 ) Estimation and adjustment of the systematic errors requires a comprehensive and efficient adjustment or calibration technique that can remove most systematic errors during the post processing of the ALS observations (Cramer, 2002). Even with the best error adjustment or calibration methods, errors still exist in the data especially in the data collected with high flight dynamics over mountainous topography and/or relatively flat te rrain with very few features (Congalton and Green, 1999). The magnitude of such errors may not be significant for many applications. However, the positional errors remain substantially large (Schenk, 2001), which introduce significant elevation errors in t errains with high relief, and applications such as change detection, geological fault detection, stream center line detection etc. can be affected by small magnitudes of errors in positions and elevations. Ensuring the accuracy of final LiDAR output requi res extensive post processing of the raw data (Slatton et al., 2007; Wehr, 1999; Luzum et al., 2005). This is achieved either by comparing and correcting collected data with ground truth, known as absolute adjustment, or by minimizing the positional and el evation discrepancies between the overlapping areas of data strips, known as relative adjustment. Due to the discrete nature of laser scanning and its monochromatic intensity, it is hard to pin point a particular target in the point cloud data with fair ce rtainty. Therefore, in terrains with low relief, where elevation values do not change significantly within a few meters along or across the flight direction, only elevation differences between the data and ground truth can be adjusted, and biases such as G PS positional bias, bore sight angular PAGE 26 26 discrepancies, scanner scale offset and range bias cannot be estimated. Bore sight angles are the angular offsets between scanner frame and IMU body frame measured at the center of IMU body frame in X, Y, and Z direct ion s The nature of laser data requires the development of adequate algorithms to photogrammetry, the laser points sample the shape of the over correspondence is practically impossible to establish under such conditions, and therefore shape based rather than traditional point based algorithms are needed. In terrains with high relief such as hilly terrains or surfaces with a significant number of break lines, the surface geometry can be compared and matched with the overlapping surfaces or with appropriate ground control information. Figure 1 3 illustrates the relative differences in conjugate point positions and elevations from the two overlappin g strips. The comparison has been made for the flat and sloping surface (with 30% slope). One potential e parameters that implies that not all the systematic errors may be recovered independently. Therefore, the question arise s : W hich errors can be recov ered and how can they be resolved ( Filin, 2003 ) ? Relative adjustment is again divided in to two categories. One is purely based on data adjustment i.e. transforming the overlapping strips using 3D rigid body motion to fit data from different swath s to each other (Fillin, 2003; Skaloud and Lichti, 2006). The necessarily adjust all the systematic errors in the data. The other type of relative adjustment involves a mathematical model relating LiDAR system parameters and PAGE 27 27 observations to the ground point coordinates, and minimization of the relative differences (3D Euclidean distances) between the overlapping strips with respect to system parameters. In systematic error adjustment through strip matching, the relative di rection of flight of the overlapping strips plays a crucial role in emphasizing the differences in conjugate features (Habib et al., 2009). For example, the positional error differences between two overlapping strips are more pronounced when flown in oppos ite directions. Figure 1 3 illustrates the cumulative effect of errors due to aircraft attitude (roll, pitch, and yaw as 0.001deg each), range bias (0.01m) and scanner scale (0.0001) on computed laser ground points. Figure 2 2 Distribution of elevation and positional discrepancies across the overlapping region between two adjacent strips on flat surface (top) and sloped surface (bottom) PAGE 28 28 The errors were simulated for two strips flown in same and opposite directions over a flat terrain at 1000 m eters of f lying height with a maximum scan angle of 30 deg rees In general, due to the systematic discrepancies in the angular observations and scanner scale, the horizontal positional errors tend to be greater in magnitude compared to elevation errors. In addition the relative differences of the systematic errors in horizontal position do not change significantly across the overlap area, whereas it changes from positive to negative in case of elevation error. This means that minimizing elevation differences betwee n the overlapping strips resolves only a few systematic errors compared to matching the features that represent horizontal positional discontinuity. By minimizing 3D Euclidean distances between the possible conjugate features the calibration parameters ca n be estimated If the strips are flown in the same direction, the positional errors due to pitch tend to be of same magnitude and direct ion in both the strips (Figure 2 3 ) and therefore, may not be completely resolved through surface matching and require ground control information to accurately geo reference the data Since it is hard to identify individual conjugate LiDAR points in the neighboring strips, it is appropriate to match the groups of points or surface elements representing certain geometric a ttributes in the overlap area. This is known as surface matching. Optimal features or surface elements are important in strip matching as their positioning and orientation greatly affect the surface matching results. The correction of laser points by mean s of a linear transformation focuses on the e ects of the systematic errors but not on their causes. This may not always be appropriate; an analysis of the similarity transformation reveals that not all error e ects PAGE 29 29 can be modeled, and thus removed, by thi s transformation. Therefore, some accuracy may be lost and algorithms that are more complicated may be needed. The existing methods mostly focus on the height inaccuracies. Existence of planimetric o set in the data will not be fully compensated for with t his adjustment and artifacts are likely to remain (Figure 1 3) Figure 2 3 Comparison of error distribution over an overlap region (60% overlap) between two adjacent strips flown in ( t op) same and ( b ottom) opposite directions: ( l eft) a cross flight profi les of the overlap region; ( c enter) t op views of the overlap area; b lack scan line is ground truth; ( r ight) r elative e rrors between the overlapping strips: a long flight position error (X e rror in black), across flight position error (Y e rror in red) and el evation error (Z e rror in blue) Another noticeable deficiency of the available adjustment procedures is the error recovery model that is usually selected In most cases, the systematic discrepancies are model ed as rigid body transformation parameters on the laser point coordinates ( Maas, 2000 ) ; they can therefore be regarded as data driven solutions. The main reason for PAGE 30 30 this approach is that most of the time users do not have access to the system parameters. However with the advances in current technology this is not a problem any longer. Hence, t his research primarily focuses on the accuracy of the ALSM data and the improvement of the calibration process by implementing a comprehensive error recovery model Definition of Accuracy for LiDAR Data There is no standard criterion for defining the accuracy of the data obtained from LiDAR ALSM samples the topography of the terrain in such a way that it is difficult to establish a correspondence between a LiDAR located ground point and the actual position of th at point on the ground. This is because a certain amount of uncertainty is associated with the geo location of every laser point due to the presence of different types of errors in the subcomponents and there is practically no redundancy in the overlappin g regions. For the same reason s it is difficult to assess horizontal accuracy. However, t he presence of certain landmarks in the scene makes it easier to compare the horizontal position from diff erent flight lines Effective accuracy assessment requires ( 1) design and implementation of unbiased sampling procedures, (2) consistent and accurate collection of sample data and (3) rigorous comparative analysis of the sample data ( Congalton et al., 1999 ) As explained earlier assessing absolute accuracy is expe nsive and time consuming. Since the accuracy of the LiDAR points varies widely over the mapped area, it is impractical to collect GPS control points for the entire survey area. At the same time, r elative accuracy can give a fairly good idea of overall accu racy. A f ew ground control points can also be incorporated to assess the absolute vertical error bound Relative accu racy refers to the accuracy with which two overlapping ALSM PAGE 31 31 swaths conform to each other. Since relative accuracy is an important criterio n to assess the quality of LiDAR data in this dissertation, only the scope of relative accuracy will be explored. Statistical Accuracy Assessment The first and very important step in a system analysis is to specify a system error budget. It specifies a th eoretical limit o n the performance a chieved by a particular system. A detailed discussion of an error budget for LiDAR can be found in Baltsavias (1999). As mentioned earlier, the theoretical performance is generally not achievable in the field. Therefore, in order to calculate the actual field accuracy, extensive GPS surveys as conducted and the LiDAR data is compared with the collected control information (Airborne 1, 2001). However, such ground surveys are very expensive and time consuming. In addition, the number of control points collected th rough these surveys range from a few tens to thousands in few cases. This number represents a tiny fraction of the number of collected LiDAR points. Since LiDAR errors are not constant within a flight line (Appendix A) especially across the flight direction a ground control information with a few control points simply cannot provide a detailed error analysis. Many times DEMs are used to assess the accuracy of the LiDAR data. There are many studies to understand the effects of the grid size and interpolation on the accuracy of the final product It has been established through these studies that rasterization of the data degrade the accuracy of the raw LiDAR data especially at the breaklines or the height variation. A relative comparison of the raw LiDAR data in the overlapping strips gives more information about the accuracy than the other methods discussed earlier. In this dissertation, the relative accuracy between the flight lines has been assessed and used to est imate the quality of LiDAR data Relative elevation accuracy PAGE 32 32 can be computed by measuring the difference between the heights of the overlapping flight lines over the sample area. Different sample areas are selected to get an estimate of varying accuracy in the survey area (equation 2.1) (2.1) The survey area can be classified into various features such as manmade features (buildings and roads), natural features (grass and trees) etc. and the accuracy can be specified for e ach category separately. Figure 2 4 Horizontal accuracy assessments Horizontal accuracy is more difficult to assess than vertical accuracy. Figure 2 4 illustrates the difficulty in measuring the horizontal or planimetric accuracy when a feature is ro tated and translated. Therefore, using a term such as maximum or minimum accuracy or accuracy with a confidence level in defining horizontal accuracy helps in understanding the quality of the data. To illustrate the concept, few samples were collected from a LiDAR data set. These samples used for horizontal accuracy assessment are building points. An approximate outline of the border s of buildings has been obtained by fitting a line to the points falling on the edges of the building. The technique used for fitting line s to the data is explained in detail in Chapter 3. The PAGE 33 33 method is called outline tracing and regularization. The statistics obtained from the given data set ha ve been tabulated below: Table 2 1 R elative accuracy of LiDAR data from different samp les with C. L (1 ) Samples Horizontal Relative Accuracy Vertical Relative Accuracy Along Flight Line Across Flight Line (m) mean (m) x (m) mean (m) y (m) Flat Roof 0. 289 0.0648 0.24 01 0.0386 0.060043 0.042457 Gable Roof 0. 3 18 0.0 462 0. 2010 0.02 643 0.01 57 0.0 7 2 41 Road 0.071521 0.0 50573 Trees (approximate assessment) 0.15930 0.180340 PAGE 34 34 CHAPTER 3 LINEAR AND PLANAR FE ATURES BASED STRIP M ATCHING Georeferencing LiDAR D ata In order to make s ense of the points collected by ALSM system, they need to be transformed into user defined reference frames. The transformation involved in direct geo referencing of the LiDAR data is a complex process. All the reference frames related to various subcompon ents such as GPS laser scanner and IMU must be defined and their spatial relationship s with each other must be computed. Generally, manufacturer s provide the spatial orientation s and translation s between these reference frames, also known as bore sight angles and le ver arms. However, due to temperature changes and platform dynamics, these parameters need calibration, which is performed after every flight mission. Reference Frames Figure 3 1 R eference f rames for ALSM s ystem ( adapted from Schenk, 2001 ) F igure 3 1 illustrates laser scanner, GPS and IMU reference frames. For most of the geodesy work, the right hand reference system is standard. Usually IMU and laser scanner frames are oriented in the same direction A calibration is performed to obtain the misalignment a long their axes known as bore sight angles. PAGE 35 35 (a) (b) Figure 3 2 (a) WGS reference ellipsoid in relation with local level reference frame (ENU), (b) vector relationship among different sensor frames in WGS reference system General error models for an airborne LiDAR System can be written as: (3.1 ) w here is the laser point on the ground and is the location of GPS phase center. is the rotation from local (navigation) reference PAGE 36 36 frame to WGS84. is the rotation from body reference frame to local ellipsoidal reference frame. is the offset vectors between laser firing point and the phase center of the GPS antenna. is the rotation between laser altimeter and the body frame or bore sight angles is the laser sc an angle or the rotation between the laser beam and the altimeter. is the laser range and is the range bias. are the random errors introduced from various sources into the measurements. ( 3. 2) After introducing constant errors such as GPS offset INS offset misalignment error scan angle error range bias and time dependent errors such as GPS drift and INS drift equation ( 3. 1) take s the form shown in equation ( 3. 2). The a bove equation can be viewed as the function of observations such as laser range, scan angles, GPS position and INS orientation, and unknowns like GPS and INS offsets and drifts, range bias, scan angle error and misalignment angles. i.e. = = The magnitude of some of the errors can be reduced by taking care of different aspects while flying and collecting the geospatial data. For example, the range error due to divergence of the beam is affected by the geometry of the target. H eight discontinuities such as the walls of houses, or cliffs and steep slopes exaggerate the PAGE 37 37 divergence error. This problem can be minimized to a certain extent by flying low and using a small laser spot Atmospheric effects c an influence the accuracy of the laser rangefinder and become significantly more critical at higher altitudes. These atmospheric effects are wavelength dependent so they can vary in magnitude depending on the particular wavelength used in the system. The c orrection for velocity change of light in the atmosphere is given as a spherical range correction. Similarly, GPS and INS exhibit some amount of systematic errors in the form of constant offsets and time varying drifts as mentioned in equation ( 3. 2). Metho dology The discrepancies in the raw LiDAR data can be adjusted by two methods. One is to geo reference the data using ground control (absolute quality control) and the other is to adjust the data in order to match the overlapping regions of adjacent swaths (relative quality control). The concept of least square adjustment to find the camera calibration parameters ( Skaloud, 2006 ) can be extended to LiDAR or laser scanning systems. However, recovery of the LiDAR calibration parameters is much more complex due to limited information. As explained earlier, t he two major problems are the non redundant determination of laser points and the unknown correspondence between laser points and the spo t they illuminate on the ground. However, in adjacent laser swaths we h ave redundant features rather than redundant points. Therefore, instead of using tie points tie features can be used. A g able roof is one example. The roof planes can act as conjugate planar features and the break line between the two (where the planes me et) can be used as a linear feature ( Figure 3 3 ) PAGE 38 38 Figure 3 3 Tie f eatures The procedure of relative quality control can be summarized in the following steps Formulation of E rror E stimation M odel The first step towards relative quality control is to find the suitable constraint features such as points, lines or planes. So far, planar surfaces have proved to be the best features for the adjustment of the raw LiDAR data. However, all the errors in the data cannot be adjusted using these features and using a combination of different types of features includ ing horizontal flat surfaces (playgrounds, flat roofs, parking lots etc.), sloped planes (pitched roof houses) and linear objects ( break lines ) generally provides the best results The LiDAR e rror estimat ion and adjustment process is based on two underlying models : Functional and Stochastic The functional model describes the mathematical relationship between the LiDAR observations and the unknown parameters, while the stochastic model describes the statis tical characteristics of the LiDAR observations. Therefore, the stochastic model depends on the choice of the functional model ( Fillin, 2003 ) In case of planar features laser points from the adjacent swaths are constrained to lie on a plane. Therefore w e can write a functional model for the laser points as, PAGE 39 39 (3 .3 ) Let be the point we wish to lie in the plane or the c oordinates of the center of the plane, and let be a nonzero normal vector to the plane. Parameters a, b and c are the direction cosines of the normal vector. The desired plane is the set of all points or such that it satisfies equation ( 3. 3). The dot product of th e normal vector and the distance between any two points on the plane is always zero. If available, the parameters of the control surfaces ( break lines and planes, in this case) should be added in equation (3 .3 ) for absolute adjustment. However, in order t o adjust the data sets relatively, another constraint equation is introduced i.e. the sum of the squares of the direction cosines should be unity. The equation can be written as: ( 3. 4) Equation ( 3. 3) can be viewed as a function of observation s and unknowns (from equation ( 3. 2)) and additional unknown parameters of the planar surfaces i.e. S imilarly, equation ( 3. 4) is also a function of the unknown parameters i.e. h(S). E quations (3 .3 ) and ( 3. 4) of the functional model are non linear Therefore, they must be linearized before adjusting them in least square sense. The linearized forms can be written as: ( 3. 5) PAGE 40 40 ( 3. 6) w here, and are the misclosures for equation s ( 3. 5) and ( 3. 6) and are the corrections. The geo referenced laser point is a function of (indirect) observations, instrument parameters and corrections. The (indirect) observations are the laser point in the sen s or frame, the GPS position and the IMU attitudes. Instrument parameters are the relative orientation angles between the ALM sensor frame and the IMU frame, and the GPS antenna eccentricit ies They are the result of (manufactu rer) lab calibration and aircraft installation procedures, respectively. Corrections are added to account for residual systematic errors in the position and orientation data and for uncertainties in the instrument parameters. These are unknown and need to be determined. The stochastic model describes the statistical properties of the observations i.e. the correlation between the LiDAR observations. In addition, many unknown parameters may be highly correlated. In that case, it may create a singularity in n ormal equation system. Sometimes pseudo measurements are introduced instead of complete unknown quantities. In our stochastic model, we can assume that some unknown errors are uncorrelated and have zero mean. The covariance matrix of the geo referenced las er point can be computed as: (3.7) where is the c ovariance matrix of the laser point in the WGS84 reference frame the c ovariance matrix of GPS position in the WGS84 ref erence frame the PAGE 41 41 c ovariance matrix of the laser point in the sensor frame the c ovariance matrix of the IMU attitudes the Jacobian matrix for the IMU attitudes and the Jacobian matrix for the laser point in the sensor frame The covariance matrices for the GPS position and the IMU attitudes are obtained from prior processing, typically from one form of least squares. Initially the correlations between the positio n and the attitudes are neglected. The covariance matrix of the laser point with respect to the sensor frame is derived from the respective sensor model. The laser point position with respect to sensor frame is a function of the observ ed scanner angle and laser range, calibration parameters from the manufacturer and corrections to account for possible changes in the calibration parameters. The covariance matrix of the laser point in the sensor frame is given by: (3.8) where and is the c ovariance matrix of the laser point in the sensor frame, the co variance matrix of the ALSM observations, the Jacobian matrix for the ALSM observations, the s tandard deviation of the laser range observation, and is the s tandard deviation of the scan angle observation. PAGE 42 42 Figure 3 4 Incidence angle of the laser beam and surface normal of the las er footprint While adjusting for the surface parameters i.e. direction cosines of planar or linear feature s e q ual weight s are usually assigned to all the points obtained from different overlapping flight lines. H owever, if we consider Figure 3 4 this s hould not be the case. The individual points of th e dataset, even if gathered with constant n avigation and range measurement accuracy, must have a large v ariation in their accuracy. The points inside the strip overlap area are likely to have better quality with measurements that are almo st perpendicular to the terrain than those where the range measurements are affected by high incidence angles. Higher incidence angle means elongated footprint and hence, higher uncertainty associated with point location. K eep ing this in mind, higher weights should be assigned to the points with better viewing geometry than those with higher incidence angles while adjusting swaths in least square sense. However, it is difficult to come up with a weight assigning technique ba sed on the footprint shape. This is because footprint shapes are difficult to e stimate. Therefore, in this dissertation, the uncertainty due to target geometry has been ignored. PAGE 43 4 3 Least S quare A djustment of LiDAR Strips After the working model is in plac e, the next step is to obtain the constraint features (linear or planar). This is explained in detail in the next chapter. Extraction of tie surfaces in the adjacent swaths is a crucial step towards correction of LiDAR data. Finally, systematic errors are recovered using the Gauss Helmert model for least square adjustment. After constructing the functional and the stochastic model s for the data, the next step is to look for conjugate features (roof planes, playgrounds, breaklines, etc.) in the overlapp ing area s of adjacent swaths. The overlap analysis identifies common areas in adjacent flight lines that meet specific criteria such as size, flatness and slope while ensuring that the overlapping areas contain a comparable number of po ints. After all of the qualifying overlap areas have been identified, the adjustment values for all the observations are computed using a least square adjustment model, which minimizes the sum of the weighted squares of the residuals. The solution can be written in the form of a set of normal equations as: ( 3. 9 ) This can be written as : (3. 10 ) and are the derivatives of equation 3.3 w i th respect to calibration parameters ( and surface parameters B is the derivative of the same equation w ith r espect t o the PAGE 44 44 observations Similarly, H is the derivative of equation 3.4 w ith respect to surface parame ters P is the weight matrix for the LiDAR observations. P c is the weight matrix for surface parameters. Equation 3.8 is solved iteratively in weighted least square sense to estimate the systematic component of LiDAR dataset. The final step in volves the assessment of the final accuracy once the corrections have been applied to the data. The residuals of the planes can verify the improvement. In order to have a more tangible idea of the improved accuracy; the adjusted data set should be compared to the ground truth. Sampling and Classification In this dissertation, t he classification of the LiDAR data has been used for two different purposes. One is to assess the relative accuracy at different sample areas and the other is to use the sampled dat a (planar surfaces and the breaklines) for systematic error adjustment. The accuracy of the LiDAR data varies widely over a terrain depending on various factors. These factors include scan angle, slope of terrain, surface roughness and reflective character istics of target, flying height, locally systematic and random errors etc. As mentioned earlier, s can angle and the slope of the target in the terrain determine the shape or the spread of the laser footprint and are major factors to be considered wh en asse ssing the accuracy of the data. Other important factor is the surface roughness Surfaces with large roughness values are likely to have lesser accuracy compared to less rough surfaces. Since manmade features are generally smoother ( with well defined geome try) than the natural objects, surface roughness can be used to classify the LiDAR data. The details of the classification process are explained in the subsequent narrations. PAGE 45 45 The full sample set is partitioned into separate sample sub sets for each class, A generic sample set will simply be denoted by D. Samples in each sub set D j are assumed to be independent and i dentically distributed (i.i.d.) according to some true probability law p( x  wi ). Therefore, the posterior probability according to Bayesian Rule is given as: ( 3 .1 1 ) The decision rule for Bayesian classifier for 2 class case would be: choose w 1 : if P ( x w 1 ) P(w 1 ) > P( xw 2 ) P(w 2 ) choose w 2 : otherwise and P( Errorx ) = min [ P( w 1 x ) P( w 2 x ) ] To design an optimal classi fier we need priors P(w i ) and likelihood density distribution p( xw i ), but usually we do not know them. The solution is to use training data to estimate the unknown probabilities. PDF Estimation with Parzen Window For distribution computation a Parzen Win dowing scheme was used. According to this, after defining a grid of some resolution, we can compute the distribution value of each d dimensional feature point (x) of the grid function. (3.12 ) The size o f the parzen window was computed by using the following equation: where n is the amount of training points for each class and h1 is an initial window size. PAGE 46 46 (3.13 ) w here n is the number of d di mensional feature points, is the volume of the parzen window, x i Gaussian Figure 3 5 PDFs of mean of first stop elevation and surface roughness of non ground points Figure 3 6 Study area for accuracy assessment: classifie d LiDAR points based on surface roughness (brown: higher values of roughness like forest; green: lower values of roughness like planar rooftops) PAGE 47 47 Figure 3 5 shows the density distribution of first stop elevation and the surface roughness for trees and build ings. The surface roughness is defined as the standard deviation of the heights of LiDAR points compared to the neighboring points. The ground has been excluded from the study area. Figure 3 6 shows the classified area for rough surfaces such as trees and planar surface such as building rooftops. Sample Collection In order to assess separately, the classified points are further extracted from the scene. The objective of this classification is to obtain as homogeneous a set of sample s as possible. The accura cy depends on the topography of ground. In order to get a more comprehensive idea, the data were classified according to the topography. The study area (Hogtown) was for its diverse topography i. e. with a lot of vegetation and manmade features such as b uildings and roads As mentioned in Chapter 2, a statistical analysis was performed on the samples to estimate the relative accuracy. The results were tabulated in Chapter 2 (Table 2 1). (a) (b) PAGE 48 48 (c) (d) Figure 3 7 Extracted planar features: (a) extracted buildings (b) gable roof (c) cross section of gable roof (d) complex gable roof ( colors indicate different flight lines) (a) (b) Figure 3 8 Extracted features: (a) a road classi fied based on elevation (b) small road patch shown in two overlapping flight lines. PAGE 49 49 Efficient Segmentation of the LiDAR Data F or error adjustment purpose, a MATLAB routine was written to classify the planar surfaces in the LiDAR data. It is based on the approach explained earlier in this c hapter. The raw LiDAR points have been classified based on surface roughness parameters i.e. the standard deviation of elevation and the rate change of surface normals The planar surfaces are less rough than the trees o r shrubs etc. Similarly, the surface normal would not change rapidly over the planar surfaces compared to rough of maximum likelihood. Figure 3 8 illustrate the results Extraction of the Build ing Points Because the performance of the calibration algorithm depends largely on the accuracy of classification, only non occluded planar and linear surface s have been pulled from the extracted data. A region growing clustering algorithm is implemented t o isolate each of the detected planar surfaces, including roads. Further, the best planar surfaces have been chosen manually for input to the adjustment model for data calibration. This step is not automated and more work is needed in the future to make th is algorithm fully automated. Figure 3 9 illustrates a few of the extracted building planes. PAGE 50 50 Figure 3 9 Extracted building with planar surfaces Boundary Extraction and Tracing Because boundaries or break lines are important features in matching the overlapping strips in the estimation process, it is important to extract these boundaries with precision. Figure 3.10 highlights the best orientation of the linear features against the scanning pattern. The laser shots th at hit the boundaries of the man made structures are desirable as boundary points (Caceres, 2008). Once the roof planes of the buildings have been extracted the next step is to find the points lying on the edges or outline of the building. This is known as boundary tracing (Sampath and Shan, 2007) In computer vision, the convex hull is a well established method to extract the shapes or boundaries. The same method is employed in this dissertation but with a little modification as explained by Sampath and Shan PAGE 51 51 (2007) Convex hull works by finding the neighboring point or pixel that make smallest angle with th e reference direction (Figure 3 11) This is a great technique for the shapes that have convex boundaries Figure 3 10 Extracted building with planar surfaces Fig ure 3 11 Modified convex hull (Sampath and Shan, 2007) PAGE 52 52 Figure 3 1 2 Traced boundaries of a building roof However, in the real world, the shape of the buildings can be both convex and non convex ( Figure 3 1 2 ) To accommoda te for the non convex shapes, a slight modification is made in the regular convex hull algorithm. While searching for immediate neighbors, a search radius is specified. In case of LiDAR points, the spacing of the points along the scan lines is smaller than across the scan lines. Therefore, instead of a search radius a search rectangle of size 30cm x 90cm is used. This modification captured the outlines of the building very well. Regularization and Least Square Adjustment This technique is similar to what Sampath and Shan (2007) has explained. Factors to be considered while implementing least square regularization include: The longer sides of the buildings mostly define the shape of the building. Therefore more weight is given to the longer edges of the bu ildings. The edges or sides of the buildings that have better incidence angle of the laser beam would have better accuracy. Hence, the incidence angles are computed for the border points and higher weight is given to the side with better geometry with the laser beam. Since the spacing of the points along the scan line is much less than across the line, the edges that are either parallel or perpendicular to the scan lines will have better PAGE 53 53 accuracy than the edges at certain angle s with the scan line. This fac tor has been considered while giving weights to the points as well. Evaluation and Results Hogtown, a residential area near the University of Florida campus was chosen to evaluate the proposed algorithm. A small planar ground was chosen as one of the cal ibration area s The second calibration area was chosen over a group of buildings with gable rooftops oriented in different directions As explained earlier, planar patches of slopes in different directions are useful in retriev ing the systematic offsets in Roll and Pitch. The data were collected with ALTM Gemini LiDAR system at 100 kHz PRF and approximately 700 m of height above ground. The maximum scan angle was 20 and the trajectory positions were captured at 10Hz with the dual frequency Trimble/Ashtec h GPS receiver. The orientation measurements of the trajectory were achieved at a rate of 200Hz using LN200 tactical grade IMU. The developed method for LiDAR calibration was implemented in a M ATLAB environment and simulated data was used to test the adjus tment model and to assess its performance. The a ircraft trajectory was processed using POSPAC, a GPS and IMU data processing software that implements a Kalman Filter to optimally combine the measurements from GPS and IMU sensors. The range measurements fro m the laser sensor were converted into ground positions in UTM coordinate system using aircraft trajectory (sbet) in Optech LiDAR data processing software, DASHMAP. The GPS conditions were almost optimal for the data set, with favorable satellite geometry and close separation from the reference receiver. The flight lines were short. R apid change s in direction (sharp turns in between the flight lines) and the increased horizontal acceleration s make the individual systematic errors within the inertial system well observable by GPS position and velocity updates while operating PAGE 54 54 flight lines of short duration limits their accumulation. These facts, together with the low flying height, contributed to the good determination of the GPS/INS trajectory (RMS of X, Y a nd Z coordinates was within 7 cm and orientation was less than 0.0 0 5 ). For the purpose of calibration, a small planar patch (approximately 15 m long and 1 0 m wide) from a road (Figure 3 14) and six roof planes with different slope orientations were selected from the classified data set Due to the high volume of point cloud, only three flight lines were used to calibrate the data set. In addition, the data set was evaluated for systematic orientation error i.e. roll, pitch and heading. Given the cross correlation among other error par ameters, their simultaneous estimation is difficult. The initial values of orientation offsets were assigned as zero. The algorithm converged in seven iterations for these initial values. Table 3 1 depicts the systematic error estimates for the given data set The standard deviation associated with each estimated parameters, points towards the reliability of the estimation process. Table 3 1 S ystematic error estimates for H ogtown dataset Estimated Error Parameter s Estimated Parameter Value s ( in degree) Sta ndard Deviation ( in degree) Roll Error 0.128 0.0009 Pitch Error 0.0503 0.0026 Heading Error 0.0013 0.0 0 16 After the calibration, the method has successfully reduced residual errors in the data set. The relative accuracy of the adjusted data is tabulated in Table 2 1 in Chapter 2 The proposed model can be extended in a straightforward manner to include other systematic errors either in LiDAR system or in navigation data; however, such expansion may result in strong correlation among some parameters The conditions for the calibration flight should be chosen as favorabl y as possible for the GPS/INS PAGE 55 55 integration and other systematic influences in LiDAR observations sh ould be estimated by independent means whenever possibl e. Another data set from the UF Campus was teste d and compared with the planar surface only adjustment model. At first, the calibration parameters were estimated using LiDAR points from only flat surfaces (abundant in the data). Later, building boundaries or breaklines were also introduced along with th e flat surface points in the estimation process. A c omparison of the correlation matrices obtained from calibration with flat planar patches ( left) and with bounded surfaces ( right) are shown in Table 3.2 Figure 3 13 LiDAR data classified based on elevat ions (University of Florida c ampus) The specifications of the LiDAR data used in the experiment are as follows: PRF=125 KHz, scan frequency= 45 Hz, Max scan angle = 20 deg, beam = Narr ow beam (0.18mrad divergence), the n umber of flight lines used for cali bration was four. The c alibration parameters obtained in both the cases were: Case (1), only flat surface elements oriented in horizontal direction : roll error= 0.008 pitch error= 0.017, yaw or heading error = 0. 00 10, range offset = 0. 0 1 0 scanner sc ale = 1. 000 02 Case (2), flat surface elements with linear features and breaklines: roll error = 0.0 09 6, pitch error= 0.0 14 3, heading error =0.006, range offset = 0.0 42 scanner scale = 1.00072. PAGE 56 56 Figure 3 14 Extracted boundaries (University of Florida c ampus ) Table 3 2 Comparison of correlation analyses of flat planar features only and bounded surfaces The comparison of the two correlation matrices shows that the correlation among the estimated parameters is very high when only flat surface points are used in the parameter estimation. The attitude errors such as roll, pitch and heading are difficult to de correlate when there is no slope change present in the sample points. Planar surfaces oriented in various directions so well in estimating these para meters. However, in my scenes, such well oriented planar surfaces are not present. In those cases, the breaklines or the boundaries of the manmade features can be utilized to achieve similar PAGE 57 57 results. The breaklines define the horizontal discontinuity and h ence are very useful in matching the overlapping strips accurately. By minimizing the horizontal discrepancies, the attitude errors can be estimated with great reliability as apparent from the cross correlation matrix in T able 3 2. From the above results, it can be summarized that a systematic error adjustment method based on feature matching especially bounded surfaces is more efficient in estimating the systematic biases in the data than the planar surfaces only method. Inclusion of boundaries along wit h planar surfaces render a more stable solution than the other proposed methods using only planar surfaces (by Skaloud and Lichti, 2006, Friess, 2006). The major drawback of the breakline based method is that the algorithm is sensitive to noise. If the fea ture boundaries are occluded (buildings under trees), the estimates will not be as reliable or stable. In such cases, these breaklines should either be avoided completely or a better breakline extraction algorithm should be used as proposed by Caceres (200 8) using spin images. PAGE 58 58 CHAPTER 4 DATA SAMPLING AND FEATURE EXTRACTI ON FOR SURFACE MATCHING As illustrated in Figure 1 2 in Chapter 1 heading and pitch errors distort the reconstructed surface as non linear warp and the magnitude of discrepancies in positi on are greater than in elevation compared to roll error and range bias. These distortions are usually along and across the flight direction. Instantaneous flight direction or local frame of reference for the system at a given moment could be different from the main flight or strip direction. As a result, each overlapping (reconstructed) surface is distorted differently. These significantly non linear distortions pose a problem in achieving good results with surface matching that rel y on a rigid body transfo rmation of one surface to another. To reliably estimate the systematic errors, the surface matching needs a comprehensive mathem atical model similar to Schenk (2003) Fillin (2003), Skaloud and Lichti (2006) and Friess (2006) Certain features that restric t the rigid body motion or six degrees of freedom i.e. translation and rotation along X, Y and Z axes (as shown in Figure 4 1( a ) ) or along and across the overlapping strips are important for surface matching. These features can be viewed as a locking syst em between the two overlapping surfaces that rigidly fits one surface to the other and restrict s them from moving relative to one another in any possible direction. The i terative closest point (ICP) method has been extensively researched as one of the sur face matching method s for overlapping point cloud data in computer vision where the size of the point cloud is much smaller and point density is more uniform than in airborne LiDAR point cloud s Sampling is one approach to handl ing large amounts of data ef ficiently in implement ing IC P as a surface matching method. Gelfand et al. (2003) and Ruckenwicz and Levoy (2001) have demonstrated that a good sampling strategy PAGE 59 59 can help achiev e the best fit between two overlapping surfaces as well as avoid local minima during convergence. Gelfand et al. (2003) explained that certain features, when chosen as samples create steep error landscape. Since ICP relies on minimizing error in least square sense, this helps the algorithm converge faster and in most cases, to the global minima. The researchers proposed a strategy to identify points that stabilize the surface matching algorithm using the covariance matrix of uniformly sampled data points. Ruckenwicz and Levoy (2001) proposed a normal space sampling method based on s electing points with variable surface normals along the entire surface. In this c hapter we propose using curvature change as the sampling criterion by selecting the points that fall at ridges and valleys on the surface. It can also be perceived as a hig h frequency filter. However, a smoothing function has also been used to suppress the noisy points that may adversely affect the surface matching process. Uniform and random sampling methods are the two commonly used sampling methods for LiDAR data. We have compared the curvature based sampling with these two methods for different terrains. The analysis has shown that for LiDAR data with relatively small systematic errors, especially scanner scale and range error, the surface warping of the regenerated LiDAR surface can be ignored and a simple surface matching based on conformal or affine transformation among the overlapping surfaces can be used to adjust these errors. For the purpose of this c hapter we have analyzed the performance of all three sampling met hods used in a simple conformal transformation based on surface matching PAGE 60 60 (a) (b) Figure 4 1 (a) Six d egrees of freedom in matching overlapping strips; (b) two ring n eighborhood Feature Selection LiD AR p oint cloud data usually contain noise introduced from system and data collection process es In surface matching, two overlapping flight strips of the same surface may not match properly due to the inconsistencies created by the noise in the data. Since the magnitude and pattern of the noise is dependent on the type of instrument and operati ng conditions, it is difficult to develop a unique feature based method to match the conjugate surfaces. However, the underlying major surface features of the topogra phy obviously are the same in the overlapping areas of adjacent strips. Our approach is to extract those features that can be used to reduce the translation and rotation instability between these overlapping strips. It is well known that the point density is not the same across the strip with LiDAR data scanners; a saw tooth scanner in this case. In the direction of fli ght, t he points are spaced uniformly near the center of the data strips than at the edges. Delaunay triangulation or polygon meshes are chos en to represent the irregularity of the point cloud data (Ghatzke and Grimm, 2006) Polygon mesh noise is another problem that destabilizes the ICP algorithm. To minimize the noise a filter similar to the median filter in image processing is implemented t o filter out the noisy points in the data. x z y PAGE 61 61 Polynomial Fitting and Curvature Estimation In a given terrain, the scale of the features and their respective curvature values vary widely. A terrain that may look smooth on a larger scale might indeed be full of subtle local features. The inclusion of features with different curvature directions ascertains the stability of the matching solution. Triangulated meshes are created for all the strips involved in the matching using Delaunay triangulation. As shown in F igure 4 1 ( b ) for each vertex point V i of the mesh, immediate neighboring points N j 1 known as first ring neighbors are identified. From first r ing neighbors their immediate neighbors N k 1 are searched for. These are known as second ring neig hbors to the ve rtex in question Using these second ring points N j 1 and N k 2 around the vertex a polynomial surface is fitted. A second order parametric polynomial or quadratic surface is fitted to these local points. The noisy points are rejected based on their elevation compared to the median of the elevations from all the fi r st ring neighbors. Inclusion of second ring neighbors will result in more smoothing. Before the computations begin for estimating the coefficients of the fit, the data associated with the selected pa tch need to be parameterized. It means that the global coordinates are represented in a local reference frame associated with each patch. The coordinates of global reference system can be represented in parameterized form as, ( 4. 1) where and The coordinates x, y and z are functions of local parameters u and v In other words, z is also a function of x(u ,v) and y(u, v). Once the local reference system is in PAGE 62 62 place with one normal a nd two tangent directions at the vertex, it is possible to write a general equation for second degree polynomial fit or bi quadratic surface fit as follows, o r ( 4. 2) The coordinates x an d y are measured along tangent direction s and z is measured along the normal direction. Term is the error term. A minimum of five neighboring points are required for the least squares solution. The cost function for the least squares fit is: or ( 4. 3) The least squares solution of equation ( 4. 1) involves finding the six coefficients of the polynomial, including constant h, by minimizing the cost function i.e. the cumulative squared error in the d ata set. After computing the coefficients of the polynomial surface, the Gaussian and Mean curvatures are computed as : ( 4. 4) w here , and As explained above, the curvatures are computed from the derivatives of the parametric surface elements. This can also be viewed as fitting a piecewise polynomial to the irregular 3D surface and computing t he curvature at each of those pieces. Then the surface is classified based on the value and the sign of the curvature. Every point is PAGE 63 63 classified into three main categories, namely valley or pit, ridge or peak, and flat surface. This procedure is repeated f or all the points in the strip. Once a point is classified, a voting procedure similar to Hough transform is followed. Each classified vertex V i casts votes for all its first neighbors N j 1 weighted by its distance from those neighbors. Each point can take three votes as valley ridge or flat. Finally when all the points are voted they can be chosen based on the highest votes received in ridge or valley or flat category. This procedure eliminates the isolated outliers or noisy points and ensures that the samples belong to a significant feature and not a superficial one The points with maximum votes in any of the categories are classified as that type of surface. The number of points sampled in a particular category can be c ontrolled depending on the numbe r of votes required by a point to be included in that category. Figures 4 2 and 4 3 illustrate the results of the sampling method based on curvature change. The number of points belonging to ridges or peaks can be down sampled by increasing the number of votes required to qualify for the category. It should be noted that the sampling process is based on the sign of the curvature and not the magnitude of the curvature. This saves some computation as well as helps in choosing the points that have small curv ature changes on seemingly flat surfaces. Surface Matching Surface matching is an important part of the data registration/geo referencing and change detection. Generally, the common features from two or more datasets are detected, extracted and matched in a least squares sense to compute the best match between the two datasets. ICP is another method to efficiently match the two datasets based on their surface characteristics (Besl, 1992). PAGE 64 64 detection and extraction of specific fea tures. As the name suggests, the closest points from two or more data sets are matched by minimizing the Euclidean distance between those points. For the purpose of this Chapter we chose the Chen and Medioni (1991) point to plane method instead of point t o point ICP. This method works by minimizing the distance between a point from one surface to the plane on the adjacent surface and is faster compared to point to point ICP. Out of two overlapping surfaces, one is used to collect the samples and the other is used as a triangular mesh surface, where each triangle is considered as a planar patch. We use the sampled points from one surface and compute their distance from conjugate triangular patches in the adjacent surface. The squared sum of all these distanc es is minimized with respect to the transformation study (1991) can be referred to for details. Results The curvature based sampling was tested on a synthetic fractal surface. Th e resulting classifica tion is shown in Figure 4 2 The points in red belong to the ridges or peaks and the points in blue belong to the valleys or pits. F eatures such as ridges and valleys are similar to break lines and can make the surface matching process more accurate by mat ching these positional variations in overlapping surfaces. There are two other types of sampling methods commonly used for LiDAR data sampling; random sampling and uniform sampling. Our purpose is to check the performance of curvature based sampling agains t random and uniform sampling for different types of terrains. For this purpose we have selected three types of terrains generated from LiDAR data. First type is a moderately varying surface created by very large sand dunes (Figure 4 3 ( a ) ). PAGE 65 65 Figure 4 2 P erformance of sampling algorithm on synthetic LiDAR data (fractal); valley points (blue) and ridge points (red). Another type of surface we selected is a relatively flat area at the bottom of a This terrain has very few small features such as a road and a few high elevation surfaces but is smaller in size. The last type is a hilly terrain created by lava flow. It has highly varied topography (Figure 4 3 ( c ) ). The idea is to find a suitable sampli ng method for each of the three types of surfaces. The results of curvature based sampling method are shown against the terrain images in Figure 4 3 The adjacent overlapping strips were matched for each surface type for each of the sampling methods. The r esulting root mean square error of the surface fit is tabulated in Table 4 1. In addition, the number of iterations for achieving the best fit for each surface is graphically shown for each sampling method. This type of analysis is useful in picking a suit able method of sampling for a given surface type. PAGE 66 66 (a) (b) PAGE 67 67 (c) Figure 4 3 Results of curvature based sampling on various terrain types : (a) sand dunes (b) flat terrain with few features and (c) mountainous terrain Table 4 1 RMSE (in meters) for different terrain types with three different sampling methods Sampling Method Hilly Terrain Sand Dunes Flat Terrain Random Sampling 1.46 0.64 0.83 Uniform Sampling 1.77 0.54 1.02 Curvature Based S ampling 0.78 0.51 0.37 PAGE 68 68 (a) (b) (c) Figure 4 4 C onvergence analysis for three different surface type s PAGE 69 69 Conclusion In this c hapter we explained that as the magnitude of the systematic errors increase (due to various factors) in LiDAR observations the non linearity of relative error di fference in the conjugate features from overlapping surfaces increases. This means that if the systematic errors in the LiDAR observations are small, a simple surface matching technique using conformal or affine transformation may be sufficient to recover e rrors in attitude and position. If the magnitude is greater, a more detailed method based on a comprehensive mathematical model is needed. In addition, it has been observed that the effect of the systematic errors is greater in position than in elevation. A sampling technique is useful for surface matching in all these cases. It reduces the computation cost and data size, helps in reaching the optimal solution and aids the automation of the surface matching algorithm. From the analysis above, it can be con cluded that the choice of a sampling method depends on the type of terrain and the accuracy desired at the expense of computer memory and computational time. For terrains with low frequency curvature variations and large features such as the surface with l arge sand dunes, any of the three sampling methods give similar results in terms of accuracy. Random or uniform sampling should be the choice for such terrains, given the fact that these two methods are simple and computationally less expensive than the cu rvature based sampling method. However, in complex terrains such as hilly surfaces, random and uniform sampling methods can easily miss crucial features necessary for a better fit during surface matching. And i n flat landscapes with minimal features, a cur vature based method is more useful to extract those important features. Therefore, the choice of a sampling method should be based on the type of terrain. PAGE 70 70 In addition, there are a few drawbacks of using ICP for surface matching as it is sensitive to data noise. In this c hapter the noise in the elevation is tackled by using a median filter however positional noise stills remains in the points. The sloping surface where the incidence angle of the laser beam is larger, positional noise in the data is immin ent. Very high positional noise makes ICP unstable and convergence to a global minimum is difficult to achieve. In such cases, a curvature based algorithm is not sufficient to achieve the best surface fit and a more detailed analysis of the data is needed. PAGE 71 71 CHAPTER 5 CALIBRATION OF LIDAR DATA USING CURVATURE BASED SURFACE MATCHING To ensure a certain level of accuracy, QA/QC is an essential part of LiDAR post processing. This is achieved either by comparing collected data with ground truth, k nown as absolute adjustment or by minimizing the positional and elevation discrepancies between the overlapping areas of data strips, known as relative adjustment. As explained in Chapter 1, d ue to the discrete nature of laser scanning and its monochromat ic intensity, it is difficult to identify the exact location of a target in the data (Huising and Pereira, 1998) This makes it necessary to formulate a strip matching method that relies on surface matching (Pfeifer et al., 2005). A comparison in position and elevation can be established if the same target with clear break lines or smooth surface is scanned by a LiDAR instrument from two different positions. Therefore, estimating systematic discrepancies through overlapping strip matching is currently a pop ular c hoice in the research community. However, large volumes of LiDAR data could be a problem when it comes to implementing a surface or strip matching algorithm. Sampling is very useful when it comes to register ing large data sets, especially airborne Li DAR The projects are usually comprise d of a few miles of long strips and some gigabytes of data are collected including GPS and INS measurements. To process and use full dataset s for adjusting the systematic d iscrepancies is a huge challenge. The alternat ive is to use only certain areas for calibration such as cross strips. The biggest disadvantage with determining calibration parameters from cross strips is that not all the errors are recoverable from such small area s of an entire project, and remaining unadjusted errors show up in other sections of the strips. As explained in Chapter 4, a nother way is to randomly sample the entire data and take a certain PAGE 72 72 percentage of the samples for point to surface matching and estimat ing calibration parameters. This way we include points throughout the entire data set and it is more likely that the retrieved parameters w ill better adjust the errors over the entire survey area. The major drawback with random and uniform sampling approach is that v aluable information is missed in terms of high frequency geometric features available in the data. The prudent approach would be to detect such features and use them to estimate the discrepancies between overlapping strips. The effect of LiDAR system errors is non linear. This means that the effect is more pronounced in certain places in a strip than in other s If such areas can be identified, their features will render more robust and optimal estimates of such systematic errors. In addition, these geometrical features can be accuracy of the LiDAR data. Accuracy is generally measured in terms of root mean square of the discrepancies between overlapping strips. ASPRS has developed guideline s for estimating the LiDAR accuracy A go od sampling strategy must have well dis tributed s amples along and across the flight direction that restrain rigid body movements in all major directions and reduce the correlation among estimation parameters. The re here along the strip is because t he calibration parameters like roll are retrieved as a transformation between overlapping strips, without involving a proper mathematical or functional model. Since er rors induced by various systematic discrepancies at different levels of the system are along the strip. In addition not all the calibration parameters are retrieved du ring the PAGE 73 73 post processing of LiDAR data. Since there is always some correlation among different parameters, applying one correction and missing other s adjustment of the overlapping strips. The r emain ing effects of such errors are more pronounced in high relief terrains. Since the method is based on minimizing 3D Euclidean distances between the overlapping strips, it is necessary to pick sample points from the places where these discrepancies are maximum. Error analysis of the overlappi ng strips flown in opposite directions helps with the same. For this study, i t is assumed that the strips are flown in opposite directions. The model can be extended to include cross strips and the strips flown in the same direction, such as the race trac k flying pattern. This method requires the calculation of the overlap extent between strips and looks for points that have the possibility of maximum spatial discrepancy. Once those areas are detected, points are picked based on their geometry to firmly l ock the strips in place This also means that the sampling method creates a steep error landscape and helps the algorithm avoid converging on local minima The other thing that helps to avoid local minima is point to plane ICP which is alternatively named as surface matching (Gruen, Acka, 2005). Introduction of points with at least 3D discrepancies between the overlapping strips or cross strips results in a low rank design matrix that further results in unreliable estimation of calibration parameters. Bas ed on the analysis done in Chapter 4, we choose curvature change based on the sampling method as our underlying surface element classification method for strip adjustment. PAGE 74 74 Figure 5 1 LiDAR data strip illustrating non uniform point density along the fli ght direction; an example of data collected with high flight dynamics. Figure 5 2 LiDAR data strip with uniform point density along the flight direction implying low flight dynamics. Proposed Method The propose d method is useful for the adjustment of ove rlapping strips where some definitive man made features are not present and data volume is very high. If the ground control information is available inside the project, the proposed method can be PAGE 75 75 used as a reliable absolute registration tool. This method r equires LiDAR data and trajectory information along with system parameters such as lever arms. The overlapping strips are matched to the neighboring strips simultaneously as done in simultaneous block adjustment in photogrammetry. The points from one strip are matched with the surface of neighboring overlapping strip s This is a modified version of to plane ICP. This method also has an advantage over the point to plane algorithm as it has a stochastic model and the measure of relia bility of the estimated parameters can be computed from the covariance matrix. The rigid body transformation parameters in point to plane ICP are replaced by unknown systematic error corrections. It is assumed that if the correct systematic corrections ar e applied to the data the strips will completely match in the overlapping regions. The algorithm uses hard surface or terrain points to find the matching surface. Vegetation or other non permanent features like vehicles in the scene can introduce errors o r biases in surface matching. Therefore, a filtering method is employed to extract ground points from the data. This method is used for adjusting relatively small survey flights such that the navigation systematic errors will not change significantly over survey time. LiDAR Equation Most recently, Triglav Cekada et al. (2009) have proposed a LiDAR equation that accounts for the effect of sensor target geometry along with the other error sources explained by Schenk (2001). For the purpose of this study, we keep our error model simple and only consider the significant errors such as errors due to position and attitude (or bore sight angles) of the platform, laser range and scanner scale that affect the LiDAR data. In this work taken into account. As PAGE 76 76 mentioned in Chapter 3, a general LiDAR equation can be written as (Skaloud and Lichti, 2006) : ( 5. 1) ( 5. 2) w here is the ground location of a laser point m and are the location and systematic error in location of the IMU center m is the rotation from local level (navigation) frame m and are the rotation and rotational offset from body is the offset vector between laser firing point and the IMU center. is the rotation between laser altimeter and IMU body frame. and are the laser scan angle (or rotation between laser beam and laser altimeter) and its systematic offset due to s canner scale. is the laser range and is the range bias. are the random errors introduced from various sources into the measurements. Surface Matching As discussed earlier, the analysis of commonly occurring systematic errors in LiDAR data reveal that the positional displacements between conjugate features in the overlapping strips are greater than elevation discrepancies. Also, the positional discrepancies are almost the same and evenly di stributed along and across the flight strip. The relative error due to bore sight roll angle between the overlapping strips is much more visible at the edges than at the center of the overlap (Figure 1 2 ). A PAGE 77 77 systematic boresight pitch error is not detectab le in the relative comparison of elevations of conjugate features or points from the overlapping strips. However, positional errors are exaggerated along the direction of the flight Similarly, errors in heading are more visible in positional comparisons o f conjugate features along the flight direction than elevation comparison. In addition, the laser points are irregularly distributed at the edges of the flight strips (in case of saw tooth scanning pattern) compared to the center of the strip. Therefore, i t is advisable to select the features other than those containing edge points. Filin (2003), Friess (2006) and Skaloud and Lichti (2006) have presented similar techniques for estimating the systematic biases in the LiDAR data that use conjugate planar s urface patches from overlapping surfaces and matches them while minimizing the error budget in the data. Generally, fitting a plane to the point cloud ensures error compensation in only one direction, i.e. the direction perpendicular to the plane. If the p lanar surfaces are oriented in such a way that the Z axis of the coordinate system coincides with the normal of the plane, the errors in X and Y axes remain uncompensated. This is usually true with flat roofs, playgrounds or road surfaces etc. For this rea son, gable roof structures are more desirable in the scene where the systematic errors are needed to be estimated and adjusted. The surfaces that lack such features, present a problem for these types of techniques. Another problem is the scanner scale erro r that tends to warp the surfaces. As a result the profile of a leveled surface might appear as a smile or frown. Analysis has revealed that the horizontal warping is greater in magnitude than the vertical error due to scanner scale. These nonlinear discre pancies (if present in significant magnitude) in the data, cannot be PAGE 78 78 compensated by simple planar surface matching or a surface matching technique without a comprehensive mathematical model. Chen and Medioni have proposed the earliest form of point to pla ne ICP. For our strips of a (ground) surface as g 1 k X The surface is defined in a three dimension al (3D) Cartesian coordinate system. i (X). An equatio n can be written for every overlap that satisfies a pair wise matching: ( 5. 3a) ( 5. 3b) In the conventional approach, one of the strips is assume d as the model surface and the adjacent strip is considered as data. A transformation is computed between the two such that it minimizes a goal function, which measures the sum of the squares of the Euclidean distances between them. In our approach, we tra nsform both the overlapping surfaces i.e. model and data by applying systematic error corrections. The ground coordinates of the laser points are expressed as functions of system p arameters, observations, systematic errors and random noise. The cost functi on i.e. sum of the squares of the point to plane distances is minimized with respect to the unknown systematic errors. The points to be matched are chosen from both the overlapping strips and almost equal in numbers. Equation ( 5. 3) represents a non linear surface that can be linearized as PAGE 79 79 ( 5. 4) where, are the estimates of surface normals, are the sample i are the centroids of the cor responding (conjugate) Since unknown parameters cannot be separated from each other, equation ( 5. 4) is non linear and can be expressed as g( l+v x) are the residuals i n the observations are the unknown systematic errors. Equation ( 5. 4) can be linearized as: ( 5. 5) or ( 5. 6) where, For each overlapping strip i and j we o btain a set of linear equations (eqn. 6) for each conjugate point surface patch pair on strip i and j A= is design or jacobian matrix x and B= is the derivative of equation (4) with respect to observations and w = x in least square sense. If ground control information is available it can be introduced as, PAGE 80 80 ( 5. 7) ( 5. 8) or ( 5. 9) where are the ground control points. If equation ( 5. 7) can be represented as equation ( 5. 8) will be the linearized form of it. Equation ( 5. 7) is the implementation of point to point ICP, where 3D Euclidean distances are minimized between conjugate point pairs instead of point to plane distances. These control points are collected either stationary or on a moving vehicle with GPS onboard over a leveled surface such as a road inside the project area. These points help register LiDAR data to the global coordinate system and allow verifying the absolute accuracy in terms of elevation. In the absence of control information the strips can be relatively adjusted among themselves using only equation (5 .5 ). The large number of points in the LiDAR data provides high redundancy to solve the system of equations ( 5. 6 & 5. 9), but t his also results in difficulty in finding the correct estimators for the transformation. Least squares principle satisfies the criterion that the unique set of computed estimators is as close as possible to the observations, considering their stochastic ch aracteristics. Least squares solution is therefore well suited for eliminating small discrepancies left in LiDAR data strips that have already been registered using INS GPS observations. The Gauss Helmert least square solutions for a system of weighted obs ervations can be given as: ( 5. 10) PAGE 81 81 where W are the weights of the observations i.e. the point coordinates (x, y, z) of the LiDAR strips. C is the covariance matrix of the observations in the LiDAR equation. As ex plained earlier, weight matrices W 1 and W 2 are the weights associated with the points (x, y, z) of the LiDAR strips. Covariance matrix C can be computed from a priori variances of the observations. These variances can usually be found in the trajectory dat a and system information provided by the manufacture of the system. Since there are 8 observation (3 for GPS, 3 for IMU, 2 for ALS observations: scan angle and range measurement), per LiDAR point, C is 8x8 matrix for one point surface patch or point point (in case of ground control information) pair. For the entire project, C will be a square block diagonal matrix of size (8*K), where K is the total n umber of conjugate point in the observations are assumed to be zero mea n and uncorrelated among them. If cross correlation information of the observations is available e.g. from GPS INS integration model, it can be used here. Similar to the Skaloud and Lichti (2006) stochastic model, we have not included the cross correlatio n terms for simplicity. ( 5. 11) where, a diagonal matrix of inverse of variances for each observation of LiDAR system. The performance of the algorithm heavily depends on the geometric character istics of input data as well as the correlation among estimation parameters. PAGE 82 82 Constrain t of the surface by features can also be seen as de correlating most of the unknown parameters. If the surface lacks certain features required to constrain its rigid body motions, the design matrix A tends to be co linear. It should be noted that even with the best geometrical features some of the parameters cannot be separated independently as they remain highly correlated. The biggest example s are INS angular bias and b ore sight angles. Also bore sight and INS roll s are correlated with the scan angle of the laser scanner. To avoid such situations, only one set of parameters is kept as unknowns in the adjustment equation. For example in above setup we keep INS angular bi ases as unknowns and not the bore sight angles. In addition, only a certain percentage of all the points belonging to some features are used as input instead of entire strip(s). The samples should be selected based on their ability to restrict the sliding of the overlapping surfaces in all possible directions. These can be referred to as geometrically stable samples (Gelfand et al 200 3 ). This also helps in reducing the amount of time it takes to achieve the optimal results through ICP. For each sampled p oint in one strip a conjugate surface patch is searched in the adjacent overlapping strip. A given surface is approximated by small surface patches. For our purpose, the surfaces are triangulated and each triangle is considered as a surface patch. The poin t from a surface is projected on adjacent triangulated overlapping surface s and the triangle containing the point is considered as its conjugate surface patch. As mentioned above, in general cases, one strip is considered as a reference strip and all the o ther strips are matched to this either simultaneously or as a pair In our algorithm, both the strips are transformed in to one object coordinate system while minimizing the goal function with respect to systematic errors. All the strips are PAGE 83 83 adjusted simul taneously through an iterative process till certain level s of accuracy are achieved. Curvature Based Sample Selection As discussed earlier, the analysis of commonly occurring systematic errors in LiDAR data reveal that the positional displacements between conjugate features in the overlapping strips are greater than elevation discrepancies. Also, the positional discrepancies are almost the same and evenly distributed along and across the flight strip. The relative error due to bore sight roll angle between the overlapping strips is much more visible at the edges than at the center of the overlap. A systematic boresight pitch error is not detectable in the relative comparison of elevations of conjugate features or points from the overlapping strips. However, positional errors are exaggerated along the direction of the flight. Figure 5 3 S ampling b ased on c urvature c hange (a) DEM of a h illy s urface; (b) a ll c lassified ridge p oints, and (c) d own sampled r idge p oints (red) Similarly, errors in heading are also visible in a positional comparison of conjugate features along the flight direction rather than elevation comparison. In addition, the laser points are sparse at the edges of the flight strips (in case of saw tooth a) (b) ( c) PAGE 84 84 scanning pattern) along a flight direction compared to across the strip. Based on Chapter 4 analysis, the LiDAR data samples were chosen based on curvature changes. Figure 5 3 illustrates the results of th e sampling method for one of the data sets in the mountainous region The number o f points belonging to ridges or peaks can be down sampled by making the curvature threshold higher to qualify for a sample point. Normal Estimation Computation of surface normals is one of the necessary steps in the algorithm. If the triangulated mesh is a ssumed to approximate a height field, which is the graph of z = f(x, y), then estimates to the first order partial derivatives and are simple to obtain from the vertex normals. The height field is implicitly represented by F(X) = F(x The gradient vector is normal to the surface, ( 5. 1 2 ) where (n 1 n 2 n 3 ) ar e the estimates of the surface normal at (x, y, z) For an implicit surface F(x, y, z) = 0 a unit length normal vector field is, ( 5. 1 3 ) where measures the first order rate of change of F normalized by the length of the gradient. For computational purposes, the cross product of any two sides of the triangular surface patch can provide the surface normal for that patch. Outlier Removal Due to the saw tooth scanning pattern of some LiDAR scanners, the point density of laser points varies along the strip more than across the flight line. This results in ill conditioned triangles near the edges. In i ll conditioned triangles one of the angles PAGE 85 85 is either very small or very large and b ecause of the near co linearity of all three vertices the estimation of normals is not stable and reliable. If included in the estimation process, these samples destabilize the estimated parameters. We have conditioned our algorithm to detect these samples by checking the condi tioning of the selected surface patch in the point surface patch pair. This eliminates most of the edge points which tend to be noisy in general. In addition, the curvature at each vertex or data point is estimated by fitting a polynomial surface to the ne ighboring points. This reduces the contribution of the noisy points in the sampling process. Computational A spects of Weighted Least Square Matching Optimum estimates of the unknown calibration parameters can be computed by minimizing the sum of the square d Euclidean distances between corresponding closest points with respect to those parameters. This becomes a non linear optimization problem and an iterative solution is achieved through least square estimation. Appropriate weights are assigned to these cl osest point pairs based on the corresponding observations. As described in equation ( 5. 10) the weights are computed as This is based on the law of propagation of errors from observations to the ground coordinates of the laser points B is the matrix of derivatives of equation ( 5. 3) with respect to the observations and P is the variance covariance matrix of the corresponding observations. The weight assigned to each point pair is decided based on the precision of their respective obse rvations. Matrix B is a block matrix and P is usually diagonal. Therefore, even with large sizes they can be easily manipulated mathematically. PAGE 86 86 Since the rotational offset parameters are generally very small for the LiDAR data strips, the sines and cosines of the rotation matrix can be approximated. The rotation matrix body frame respectively, as: ( 5. 1 4 ) where c and s are the c osines and sines. We approximate all the small angle s such that and The rotation matrix will be ( 5. 1 5 ) For further simplification, the error in scanner rotation is assumed to be only due to scale factor. As explained in equation ( 5. 1 5 ), the error in scan angle can be mainly attributed to scanner scale ( sc ) and scan angle offset ( ). However, a very small offset between IMU and laser scanner position and orientation makes it di fficult to independently estimate scan angle offset and the roll angle offset Therefore, for simplification, we assume that = i.e. the scan angle systematic error is due to scanner scale. ( 5. 1 6 a) ( 5. 1 6 b) PAGE 87 87 Similarly, systematic error in boresight misalignment and attitude offset are highly correlated. Since they have a similar effect on the ground coordinates of laser points either one of them is assumed unknown. With the given assumptions, it simplifi es our problem in estimating attitude offset, scanner scale, and range offset. Since the method is designed for small survey area adjustment, we assume that the systematic errors remains the same over the entire survey area i.e. they do not change over tim e and therefore, the unknowns INS angular bias, scanner scale error and range offset are assumed constant throughout the survey. A work flow diagram explains the sequence of the proposed algorithm steps in detail (Figure 5 4). Figure 5 4 Work flow of the proposed surface matching and systematic error adjustment algorithm PAGE 88 88 Convergence Achieving meaningful convergence of ICP variants is a well known problem. The present study deals with the problem in different ways. Our sampling method picks the points aro und certain features that tend to tie overlapping strips rigidly. This creates a dense error landscape around these features and the algorithm converges to global minima as opposed to local minima. In addition, airborne LiDAR data are pre aligned using GPS and INS observations and the systematic errors tend to be very small in for the convergence of the algorithm. Avoiding the outliers, as explained above makes algor ithm divergence proof. However, the correlation s among unknown parameters creates a problem for a reliable estimate of the parameters. All correlated parameters cannot be introduced as unknowns simultaneously in the estimation process. Only one set should be selected as unknown. Results and Discussion Same as in Chapter 4, t his method was tested on three LiDAR data set s (Figure 5 4) with ground control information available in one of the project s The first data set contained LiDAR points collected over vo lcanic surface s in Hawaii with an Optech Gemini system at 100 kHz and at a flying height of approximately 700 800m above local ground level. This data set was particularly of interest as it was flown with high flight dynamics and thick clouds scattered ove r the survey area. Large elevation and position discrepancies were expected in t he data. The set was first processed with commercial software that adjusted for some of the systematic errors. However it still showed significant mismatches over the overlapp ing strip areas. Three overlapping strips were PAGE 89 89 chosen for estimating calibration parameters for each flight. The swath overlap was approximately 200m x 500m in size. The second data set was collected over gypsum sand dunes in New Mexico at 70 kHz with max imum scan angle of 3 0 degree s Th is data set w as collected from a relatively high altitude ( approximately 2 5 00 meters above ground ) as the sand was highly reflective and could saturate the laser sensor. The flight was less dynamic and three strips were cho sen for the adjustment and calibration parameter e s timation The third data set was collected over a meteor crater surface with an Optech Gemini system operating at 70 kHz PRF and a flying height of approximately 2000m above local ground level with a scan frequency of 40Hz and maximum scan angle of 21 degree. All three data set s were first processed using initial calibration parameters such as bore sight angles and scanner scale s from immediately previous project s In the case of the meteor crater data set t hree overlapping strips were chosen for simultaneous adjustment and estimation of systematic errors in bore sight angles, scanner scale, range bias and relative positional offset. The strip overlap was approximately 50 60% in each data set Meteor crate r and Sand Dune data strips were flown in the same direction compared to lava surface data strips that were flown in opposite directions As discussed earlier, pitch bias is difficult to resolve in such cases, and some of the ground control points were use d in the adjustment procedure as explained in equation ( 5. 9) and ( 5. 10). Ground control collected at the crater site c ould help eliminate a bias in the elevation due to systematic pitch error. However, the ground control information was not used in the adj ustment model as explained in equation 5.7. A quick analysis revealed that the road used for ground control data collection was fairly PAGE 90 90 wide and flat. Introducing ground truth in the adjustment model and matching with the LiDAR data set introduced significa nt uncertainty in the position of LiDAR points and hence the variances of the estimated parameters were very high. To avoid unreliable estimates of the parameters, the ground control information was used only to adjust the elevation discrepancy in the data Figure 5 5 Three different types of terrain: (left) mountainous lava surface; (center) sand dunes; and (right) flat surface with few features. All the data sets were filtered at first, using commercial software, to remove any unwanted points suc h as vegetation, multipath or noisy points. These points can introduce bias es in the data and hence the final estimated parameters. As explained in the previous chapter curvature based classification was performed on all the strips and the sample points w ere selected from both the strips in a given overlap. An example of the sampled point from the given data is presented in figure 5 3 (b). The site is almost flat with few features (1 4m in elevation ) present in the scene. The proposed sampling method selec ts most of the sample points near the features of interest. Over the flat a few points associated with subtle local changes. It is advantageous to have a few points from the flat area s as the point to plane algorithm freely moves over flat surface s take long to converge By having more points around constraining the feature s, we PAGE 91 91 ensure a dense error landscape that enables the algorithm t o converge to a global minimum (a) (b) Fig ure 5 6 Study area 1 : (a) intensity image (b) flight swaths After removing the outliers, the final sampled points were approximately 1 5 20 % of the total data points in all three data sets Searching for conjugate surface patches was the most time cons uming part of the algorithm. The data points in each strip were also was employed to search for the conjugate surface patches of each sampled point in its triangulated neighb oring data strip. Simultaneous adjustment enables all the data to match strips to match with each other at the same time and does not require a reference strip. The initial approximations of the unknown systematic errors were set to zero. The results for m eteor crater site were obtained after 9 iterations compared to lava surface site (12 iterations). The meteor crater site results were compared with ground truth elevation. The ground truth data are collected with a GPS antenna installed on the rooftop of a vehicle. The antenna height is measured from the base of the wheels. The vehicle is driven in the middle of a road inside the survey site. Since ground truth points cannot be identified in the LiDAR data with full confidence, a nearest neighbor approach PAGE 92 92 i s used to establish the correspondence between the two. With the given configuration, the LiDAR data points will have maximum spacing of 0.9m near the nadir and 1.7m near the edges. There is a good chance of finding a laser point for every ground control p oint within a 1m radius. Figure 5 7 Study area 2; (a) DEM of a strip; (b) LiDAR points were selected using curvature based sampling from a small area (red box on DEM ); (c) overlapping strips with ground control (line across the strips) In comparing th e results, the assumption is made that the road is flat and the within a small distance of one meter Elevation differences were c omputed by subtracting the ground truth elevation s from the nearest laser point elevation s for eac h of the three strips. Figure 5 10 shows the plots of histograms before and after the adjustment was applied to the LiDAR points. A careful analysis of the histograms of strips before adjustment reveals that the PDF resembles a rectangular shape. The spread of the elevation error (compared to ground truth) in each PAGE 93 93 strip is wide and all the values have almost equal frequency. This is only possible when a strip is inclined to the straight ground truth line. As the ground truth is collected across th e strip which is almost perpendicular to the flight direction, the histogram points towards the presence of roll error with or without scanner scale error in the data. In addition the histograms are no t centered at zero value of the elevation differences. This clearly indicates a bias in the elevation of ALS data before adjustment. After the adjustment, the PDF is almost bell shaped centered at or near zero value, indicating presence of mostly Gaussian noise and no roll error. Since there was no g round t ruth available for lava surface site the results were compared with the commercial software results for the same data strips. Figure 5 8 and 5 9 illustrate the comparison of the results. T he strips were gridded into 1m x 1m DEMs and the height differences were computed between the overlapping regions of each strip with another. Evident from figures 5 8 and 5 9 the proposed algorithm based on selective sampling with stable point pairs reduced the average elevation error significantly in both the overlaps. The results show that by using an efficient sampling strategy and adjusting the flight swaths independently it is possible to enhance the calibration results and reduce overall errors in the data. Table 5 1 Calibration parameters for lava surface data Cross Correlation Matrix Estimated Parameters Parameter Value Std. Deviation 1.00 0. 0 005 0. 22 0.0 01 0. 10 0.00 19 0.00 1 7 0. 0 05 1.00 0.0 14 0.00 2 0.0 16 0.00 07 0.00 12 0. 22 0.0 14 1.00 0.00 4 0.00 9 0.00 3 0.00 2 5 0. 0 01 0. 0 02 0.00 4 1.00 0.00 1 0 .000 5 0.000 3 0. 10 0.0 16 0.0 09 0.00 1 1.00 0.0 48 0.00 2 2 PAGE 94 94 Figure 5 8 E levation discrepancies over first overlap area. (left) estimated with co mmercial software (average elevation error = 0.060m) ; (right) with proposed calibration algorithm (average elevati on error = 0.0012m) Fig ure 5 9 Elevation discrepancies over second overlap area. (left) estimated with co mmercial software (average elevation error = 0.11m) ; (right) with proposed calibration algorithm ( average elevation error = 0.0026m) Similar results were obtained for sand dune surface in approximately 5 iterations. The results of the meteor crater site are illustrated in Table 5 2 and Figure 5 10. It should be noted that the performance of the algorit hm is better compared to the commercial software in case of data that were collected with high flight dynamics such as lava surface site. PAGE 95 95 Figure 5 10 Calibration results for three overlapping flight lines : (a), (b) and (c) are the histo grams of elevation differences between three strips and the ground truth before adjustment. (d), (e) and (f) are the histograms of elevation differences of same three strips compared to ground truth after systematic error corrections were applied to them. Table 5 2 Correlation matrix of the estimated p arameters (left) and parameter values with their respective standard deviations (right) ( m eteor crater data set) Cross Correlation Matrix Estimated Parameters Parameter Value Std. Deviation 1.00 0.16 0.09 0.02 0.001 0.00 7 0.0037 0.16 1.00 0.07 0.0006 0.003 0.006 0.0024 0.09 0.07 1.00 0.001 0.00 0.001 3 0.0008 0.02 0.0006 0.001 1.00 0.007 0 .00024 0.0 0013 0.001 0.003 0.00 0.007 1.00 0.026 0.005 (a) (b) (c) (d) (e) (f) PAGE 96 96 For the other two data sets, the results of the proposed algorithm were similar to the commercial software. Table 5 2 illustrates a good solution with low correlation among the estimated para meters e xcept for the correlation between roll and pitch error This correlation is expected because the flight lines were flown in same direction and the effect s of pitch error tend to be exactly the same in position and elevation in terms of magnitude and direction for the adjacent strips. Therefore relative matching of the adjusted using ground control to some extent as done in the experiment. The correlation be tween roll and pitch also explains the slightly high standard deviation in the estimated parameters. An overall improvement in the accuracy is evident from the results. The proposed algorithm was further analyzed for its stability in estimating the system atic errors of vari ous magnitudes. Three copies of the Meteor Crater data were created with known systematic errors intentionally added to them The estimated parameters for each data set were compared against the intentionally added systematic biases in T able 5 3 As the magnitude of the added biases increased, the standard deviation of each of the parameter increased too. In the third data set, the magnitude of the added attitude errors is high. As a result the variance of the estimated parameters is hig her than the instrument specified variance. This means that the algorithm is not as stable as it was for the other two data sets with comparatively smaller systematic biases. As the systematic biases increase the conjugate features are moved farther away in the overlapping regions making it difficult for an ICP based algorithm to converge quickly and to a global minimum. The performance of the purposed algorithm was compared with the commercial software as illustrated in Table 5 3. PAGE 97 97 Table 5 3 Stability ana lysis of proposed method by comparing intentionally added (black) and recovered systematic errors [curvature based method in ( red ); commercial software in ( blue )] (Study Area 2) Data Set Parameter Values True (Added) Parameters 0.005 0.005 0.005 0.0001 0.01 Estimated Parameters 0.00443 0.00398 0.0056 0.00014 0.001 1 Standard 0.002 2 0.0017 0.0015 0.00009 0.0071 Estimated Parameters 0. 00 51 0.00 4 76 0.0056 0.0001 1 Standard 0.00 09 0.0 0 0 56 0.0015 0.0000 4 True (Added) Parameters 0.01 0.01 0.01 0.0001 5 0.0 5 Estimated Parameters 0.0076 0.0 14 0.0092 0.00022 0.079 2 Standard 0.0058 0.0024 0.0039 0 .00067 0.045 Estimated Parameters 0.00 93 0.00 81 0.0092 0.000 1 3 Standard 0.00 16 0.002 9 0.0039 0.000 0 6 1 True (Added) Parameters 0.0 5 0.0 5 0.0 5 0.000 5 0. 10 Estimated Parameters 0.0208 0.0377 0.0165 0.0001 0 0.0003 3 Standard De 0.037 0.0165 0.046 0.00 0 9 35 0.0594 Estimated Parameters Standard PAGE 98 98 It should be noted that the performance of the proposed method is comparable to the commercial software in terms of stability. When the m agnitude of the systematic is large, the commercial software stops the estimation processes without reaching any solution (Data set #3 in Table 5 3). In addition the commercial software does not have an option to estimate the error in range. Generally, in most of the LiDAR data sets the systematic biases are within a range that an ICP algorithm can handle. However, if large magnitude errors are suspected, an initial alignment is required to bring the overlapping surfaces close enough for an ICP based method to work well. Usually, the initial estimates for the adjustment model come from the systematic error parameters that were estimated for the previous flight mission. A convergence analysis of the proposed method was also performed by testing it on three di fferent terrains. These were chosen based on their surface characteristics (Figure 5 3 ). The f irst terrain contains lava flow over a mountainous region and resembles a fractal surface. The s econd data set contains sand dunes and resembles a wave like scene The third type of terrain was our test area which is flat with few features. The root mean square error was computed for each iteration and the results were compared for all three surfaces. The plot of RMSE versus the number of iterations illustrates the performance of the algorithm for different types of geometry and hence the constraining conditions present in a particular data set (Figure 5 1 1 ). It took approximately 15 iterations for the solution to converge in mountainous terrain compared to other te rrains. The amount of noise in the data definitely affects the performance of the algorithm which is particularly noticeable in the mountainous terrain. On the other hand, the sand dune data had few systematic errors to begin with, and it therefore took only 6 iterations PAGE 99 99 for the solution to converge. The results of our test area remain the best and can be explained by the inclusion of control information in the adjustment model as opposed to other two surfaces. It should be noted that the strips in the me teor crater and sand dunes were flown in the same direction compared to lava surface where the adjacent strips were flown in opposite directions. Figure 5 1 1 Convergence analysis of proposed method over three different types of terrains Conclusion In t his c hapter a surface matching technique has been presented to estimate the systematic errors present in ALS data using surface matching techniques. It is based on identifying suitable features in a data strip and matching them to overlapping surfaces fro m the adjacent strips. All the strips are simultaneously matched and transformed while minimizing an error metric as a function of systematic errors. The problem associated with the convergence of the modified ICP algorithm is addressed in two ways. One is by selecting points from features with significant curvature change and second is by choosing appropriate unknown parameters for estimation. A sampling PAGE 100 100 method based on the geometry of the features is generally sufficient to constrain the movement between overlapping strips. Even distribution of the constraining features over the survey area can ensure the stability of the surface matching. We have successfully demonstrated the use of ICP for the simultaneous adjustment of ALS data along with the elevation adjustment using ground control information. As explained earlier, it must be noted that such methods based on minimizing mismatch between the adjacent overlapping strips can deliver better results when the strips were flown in opposite directions. Because the research involves sampling large data sets and searching for nearest neighbors, this method (based on ICP) can be very time consuming and results in a high use of computer memor y However, an optimally written program can handle the computational load involved with the method. Point search is the most time consuming part involved in the process. For our purposes, all the data are not loaded in to memory at once. Data strips are only loaded when required so that the algorithm does n o t run out of memory d uring the processing. Once the pre processing is performed for removing vegetation and unwanted points, human intervention is not required until the results are comp uted. If the sample selection can be refined further to be adaptive for any terrain type wi th an optimal number of sample points, the strip adjustment or system calibration process can be fully automated. More work is required to fully explore this possibility. Finally, if applied efficiently, a modified ICP can be used successfully to adjust Li DAR strips even when no man made features are present. PAGE 101 101 CHAPTER 6 DISSERTATION CONCLUSION AND FUTUR E WORK As required for the various applications, the research community has been involved in developing standard calibration procedure s for airborne LiDAR d ata as a post processing step. Most of the commercially availa ble software empirically adjust s the adjacent LiDAR swath s or strips, which are less rigorous and may not solve all the errors. From the exiting studies (Schenk, 2003) and the analysis of simul ated LiDAR data in this dissertation, it has been concluded that the horizontal accuracy is as vital as the vertical accuracy; especially in the irregular terrains. S econdary elevation errors result from the planimetric discrepanc ies in ALS data. Keeping t his in mind two detailed algorithms have been proposed that use the features defining the horizontal extent such as breaklines and curvature changes for matching the overlapping LiDAR strips. The first method is useful for urban areas as it uses manmade f eatures such as buildings for strip matching. This method is best suit ed for areas with g ood geometrical location and orientation of the planar features in the overlapping areas. The second method based on matching the curvature changes in the topography is better suited for natural terrains that lack manmade features. The proposed method s w ere implemented, tested, and evaluated for estimating calibration parameters such as attitude errors scanner scale and range bias Although the algorithm is evaluated only for systematic offsets in the trajectory orientation scanner scale and range parameters, it holds the potential to solve for a time dependent component of attitude error parameters. In f uture work this algorithm c ould be expanded to verify its use i n various other application s, such as surface change detection. Generally, the achievable accuracy depends on the density and the shape of PAGE 102 102 the extracted features In other words, it depends on how well surface feature surfaces such as linear, planar or cu rvature features can be extracted from the data set. Based on the analysis, t his work can be further carried on to explore the following possibilities: The performance of a curvature based adjustment method can be negatively affected by the presence of ex cessive noise in the data. Use of s pin images to extract the curvature from a noisy data set and estimating the systematic errors would be a powerful addition to this method. The procedure of extracting planar surfaces using spin images has been successful ly demonstrated by Caceres (2008) in his doctoral dissertation. A similar procedure can be implemented to de noise and extract vital manmade or natural features from the terrain data. The c urvature based method can be extended to be useful in change detect ion applications. As the analysis in Chapter 4 demonstrate s the curvature based method is useful in the surfaces where the features are scarce or in abundance. Change detection methods for such surfaces can be benefitted from this sampling approach. Furth er investigation of LiDAR is needed to explore the presence of time dependent systematic errors especially GPS and INS errors. The proposed surface matching technique for the natural terrain can also be useful in LiDAR aided navigation. LiDAR intensity c an provide useful information related to position of various features. Since the proposed method uses positional information along with the point elevations to match the overlapping strips, intensity (where available) can be an additional feature for estim ating the systematic discrepancies in the data. PAGE 103 103 APPENDIX SOURCES OF ERRORS AN D THEIR EFFECT ON LIDAR DATA The positional and altitude errors in the ALSM data can be broadly categorized into three types. Errors due to system components Errors due t o measurements Errors due to target reflectance and geometry Since ALSM is a complex assembly of different subsystems, each one of them contributes towards the gross error budget. These errors can be either due to the mechanical or electronic defects of t he system or the errors in the measurement of the observables. It mainly includes GPS bias, IMU errors, atmospheric delay, bore sight misalignment, scanner related errors and laser beam misalignments, time synchronization error etc. In addition, the topogr aphy of the terrain plays crucial role in the final accuracy of the LiDAR data. Beam Divergence One of the intrinsic properties of a LiDAR scanner that strongly influences both the point cloud resolution and the positional uncertainty is the laser beam foo tprint. The apparent location of the range observation is along the centerline of the emitted beam. However, the actual point location cannot be predicted since it could lie anywhere within the projected beam footprint. Since a ll photons travel in the same direction, the laser light is much collimated and usually low in divergence. However, at highly sloped surfaces where the beam incidence angle is very high, it results in positional uncertainty as explained earlier. The effect of divergence of the beam on the accuracy of the laser ranges depends on the flying height as the height discontinuity exaggerates the error. Two factors that can be controlled to a certain extent are flying height and the beam PAGE 104 104 divergence (some ALSM systems come with an option for a narrow or wide beam). Low flying height and a narrow laser beam will result in lesser uncertainty in the computed ground coordinates. Target Geometry and Reflectance The spread of the laser footprint is also affected by the geometry of the target F or exa mple, larger slope values of terrain increase the spread of the footprint. Flat terrain and large scan angles have a similar effect (Figure A 2 ) However, the target geometry i.e. the topography of terrain cannot be controlled. The errors due to the foot print spread are random in nature and not included in the systematic adjustment models Figure A 1 Transmitted and received laser pulse (CFD) The reflectance properties of the target affect the strength of the reflected signal from the target. Constant f ractio n discriminator (CFD) (Figure A 1 ) is designed to produce accurate timing information for the transmitted and received signals but fa ils to record the correct time for exceptionally low or high intensity values of the reflected signal. PAGE 105 105 Figure A 2 Laser foot print on ground This leads to errors in ranges and hence in the positions of ground points. Platform P osition and A ttitude Apart from the target reflectivity properties and laser beam incidence angle, the main limiting factor s are the accuracy of the platform position and the orientation derived from the carrier phase differential GPS/INS data and uncompensated effects in system calibration. GPS Positioning Use of the Global Positioning System (GPS) is an important part of LiD AR mapping. Airborne GPS systems are used in LiDAR to provide positioning information regarding the trajectory of the sensor. GPS is not an autonomous system, which means that it needs external devices like satellite constellation to compute its position, and therefore, many factors contribute to the accuracy of GPS derived positions. These sources of error include satellite geometry (PDOP), the number of satellites, orbital biases, multipath, antenna phase center modeling, integer resolution, and atmosphe ric errors. The errors in GPS directly affect the horizontal as well as vertical accuracy of the LiDAR data. PAGE 106 106 (c) Figure A 3 GPS positioning (a) absolute (b) differential (c) IMU gyros and accel erometers mounted orthogonally (Colomina, 2000) Exterior Orie ntation from Inertial Measurement Unit Knowing the correct orientation of the sensor in space is a necessary but not sufficient condition for calculating an accurate transformation from the local sensor reference frame to the Earth centered reference frame (WGS84). In practice, the orientation of the platform is recorded by an on board inertial measurement unit (IMU) that is hard mounted to the LiDAR sensor. An IMU comprises three accelerometers and three gyros arranged in an orthogonal triad (Figure A 3 ) Accelerometers measure the local gravity and acceleration vector, and gyros measure the angular rate vector experienced by IMU. Typical sampling rate of an IMU is 100 500 Hz. Post processing GPS IMU software converts measured accelerations and angular rates along with GPS positioning to position, velocity, roll, pitch and heading. Accurate measurements of the roll, pitch, and heading of the platform are required to correctly determine the pointing direction for each laser pulse. Small angular errors in roll, pitch, and yaw result in positional and height discrepancies in the LiDAR data. During the large maneuvers (usually at the beginning or the end of flight strips), the accuracy of the IMU degrades and hence degrades the accuracy of the LiDAR data. PAGE 107 107 B ore Sight Angles In most system installations, the lever arms between LiDAR IMU GPS sensors can be determined separately by independent means, although this represents certain difficulties related to the realization of the IMU body frame. On the other hand the determination of the bore sight angles is only possible in flight once the GPS/INS derived orientation becomes sufficiently accurate. These angles are usually determined by the manufacture of the ALSM system and assumed constant over time. However, t he changes in temperature and dynamics on board the airplane can change these values slightly and might need calibration over a long time span. Scan Angle The current ALSM systems use three different types of scanner s : Fiber Scanner (TopoSys) Rotating Poly gon (IGI) Oscillating Mirror (Optech and Leica) Oscillating mirror based scanners have been proven to be most efficient so far. A small angular misalignment in the scan angles translates to significant positional as well as elevation error. The index erro the Z axis of the scanner ( Schenk, 2001 ) The scan angle error is not constant; rather it changes with the scan angle. In addition, the oscillating mirror slows down at the end of the scan to rev erse the direction. This results in the accumulation of laser points at the edges of the swaths. Non linearity of the scan angle also contributes to the error budget. Laser Range A discrepancy in range translates into errors in ground coordinates of the b eam. There may be error s in the laser range measured due to time measurement error s PAGE 108 108 wrong atmospheric correction and ambiguities in the target surface, which results in range walk. Atmosphere Atmospheric effects can affect the accuracy of the laser ran gefinder and become s significantly more critical at higher altitudes. Atmospheric delay of an electromagnetic wave can be classified in to two categories: dispersive and non dispersive. The n on dispersive portion is usually associated with the troposphere w hile the dispersive potion is associated with the ionosphere. Troposphere is the lower part of the atmosphere above the 18km. The troposphere experiences temperature, pressure and humidity changes due to the weather. Sinc e these variables affect the air density and hence the speed of light, changes in tropospheric conditions result in inaccurate range measurements of the electromagnetic waves. The ionosphere is the layer of atmosphere above 50 km that consists of ionized a ir. Changes in the level of ionization affect the refractive indices of the various ionospheric layers and thus affect the travel time of the electromagnetic waves. These atmospheric effects are wavelength dependent so they can vary in magnitude depending on the particular wavelength used in the system. To accurately measure the ranges, it is necessary to correct for the atmospheric delay. Atmospheric delay depends on the refractive index along the path that the laser pulse takes through the atmosphere. Fo r the purpose of airborne laser ranging systems that do not fly at very high altitudes (less than 2000m), the t ropospheric delay does not play a significant role. The refractive index of air is a function of density. ALTM Gemini system operates at a wavele temperature of 15 degree Celsius, atmospheric pressure of 101.325 kPa and 50% of PAGE 109 109 humidity the refractive index of air comes to be 1.00027371. The correction for the refraction and velo city change of light in the atmosphere is given as a spherical range correction and needs to be calculated at both the platform location and the laser footprint ( Airborne 1, 2001 ) The one way correction to the laser ranges due to the atmospheric delay can be written as, (A 1 ) where n(z) is the refractive index of the atmosphere along the ray path, the first integral is the curved path followed by the laser pulse from the aircraft to the ground, and the second integral is the straig ht line path from the aircraft to the ground. For the zenith path, equation 2.1 can be written as, (A .2) w here z 1 is the elevation of ground point and z is the aircraft elevation. Since temperature, pressure and the humidity changes with elevation, the refractive index n(z) is the function of elevation. However, the operating elevation of UF ALTM Gemini system is usually less than or equal to 1000m ; these weather conditions do not change significantly (Fig ure A 4 ) Hence, the refracti ve index of air can be assumed constant. If the laser range calculations are done based on the speed of light in a vacuum, the accuracy of the ranges can be affected by a constant offset of 27 cm at an altitude of 1000m. PAGE 110 110 Figure A 4 Profiles of temperatur e and pressure with altitude Manufacturers of the LiDAR systems and processing software take into account the refractive index of the air and hence, such systematic error is avoided. The possible error that can be introduced is due the significant changes in the refractive index of air due to turbulent weather conditions or when the system is operated at very high altitudes. Interpolation and Coordinate System Conversion The LiDAR system can scan the terrain at PRF (pulse repetition frequency) ranging from 33 to 167 KHz. However, GPS INS derived trajectory usually has a frequency equivalent to that of an IMU (200Hz for ALTM Gemini). In order to compute the position of each laser footprint on the ground, the airplane trajectory is further interpolated to mat ch the frequency of the collected laser range data. This introduces the interpolation errors in the post processed LiDAR data. In addition, the measurements from different sensors are converted into a common reference system (generally WGS 84). Data are f urther transformed into local co ordinate systems (e.g. PAGE 111 111 UTM). Inaccurate transformation parameters affect the accuracy of the resulting datasets. However, these are not considered as major sources of error and therefore, have not been explored in detail in the existing literature. Effects of Major Systematic Errors Attitude Errors Roll e rror As per the geometry in Figure A 5 roll error would translate into horizontal error dy and vertical error dz. In Figure A 5 the The horizontal displacement of the point can be derived as: (A 3 ) (a) (b) Figure A 5 Roll ( ) error (a) profile of two overlapping strips showing displacement of a point due to roll error resulting vertical (dz) and horizontal (dy) error Similarly, the vertical displacement of the ground point can computed as: (A. 4 ) PAGE 112 112 In the absence of yaw or heading, the roll error would result in a vertical error along vertical axis Z and the horizontal error along the Y axis of the body frame (right handed IMU reference frame). The rotations follow the r ight hand thumb rule and for the correct direction of the displacements, the sign of the rotations must be included in the above equations. Pitch e rror Similar to roll, pitch error displaces the laser points on the ground in both horizontal and vertical di rections. The shift of points due to roll error is usually visible across the overlapping flight line. However, the displacement of points due to pitch error is evident along the overlapping flight lines. In Figure A 6 (a), the LiDAR system with pitch angle LiDAR system with scan horizontal error dx and vertical error dz. In Figure A 6 (b) LiDAR point in the direction of flight. (a) (b) Figure A 6 e rror (a) side view of two overlapping strips showing displacement of a point due to pitch error and (b) resulting vertical (dz) a nd horizontal (dx) error PAGE 113 113 Therefore, in the real data, pitch error is evident from irregular spacing of the scan lines in the flight direction. The horizontal displacement of the point can be derived as: ( A 5 ) S imilarly, the vertical displacement of the ground point can be computed as: ( A 6 ) In the absence of yaw or heading, the pitch error would result in vertical error along vertical axis Z and the horizontal error along the X axis of the body frame (right handed IMU reference frame). Heading e rror As seen in Figure A 7 the heading or yaw error results only in horizontal displacement of the laser point on the ground and no vertical displacement. Therefore, an across strips profile analysis cannot reveal this type of error The horizontal displacement of the point along the flight line can be derived as: ( A 7 ) Similarly, the displacement of the laser point on the ground ac ross the flight line can computed as: ( A 8 ) The a bove equations are valid in the absence of pitch. The position of the ground laser point is further changed when pitch is included in the computations. Larger va lues of roll and scan angle further exaggerate the yaw error as they affect the range The discrepanci es due to roll errors (Figure A 11 ) are clearly visible on the inclined gable roof planes. The errors due to pitch are also apparent in the scene. PAGE 114 114 (a) (b) Figure A 7 Heading or pping strips showing displacement of a point due to heading error top view of the flight line showing horizontal (dx and dy) errors However, as explained above, the errors in yaw or heading cannot be discerned from the profile, ( A 9 ) ( A 10 ) ( A 11 ) ( A 12 ) The position of the laser point on the ground in the flight line projection frame can be computed from the range scan angle and all thre e components of platform attitude as well as position PAGE 115 115 ( A .1 3 ) Scan Angle Errors The index error in the scanner has a similar effect as the roll error in exterior orientation. Non linear sc anning speed of the scanner accumulates LiDAR point s on the edges (Figure A 5) Figure A 8 Scanner mirror error Misalignment in the laser beam and the vertical axis of the laser system leads to planimetric and elevation errors in the dat a similar to bore sight misalignment and attitude errors. Figure A 9 Scanner s cale ( s c ) e rror PAGE 116 116 Range Error This type of error depends only on the scan angle i.e. the effect of range error can be seen more on the edges than anywhere else. It affects b oth planimetric and elevation coordinates. Figure A 10 Range b ias ( ) Figure A 11 A cross flight line profile of two overlapping strips over a gable roof PAGE 117 117 Figure A 1 2 Effects of systematic attitude biases (roll, pitch and heading) on the n ominal surface (blue) are visible in the modified surface (red) PAGE 118 118 Figure A 1 3 Effects of systematic scan angle offset and scale bias on the nominal surface (blue) are visible in the modified surface (red) Figure A 1 4 Effects of systematic range bias on the nominal surface (blue) are visible in the modified surface (red) PAGE 119 119 LIST OF REFERENCES Akca, D., Gruen, A. (2005). Fast correspondence search for 3D surface matching. International Archives of the Photogrammetry, Remote Sensing and Spatial Informat ion S ciences 36 (3/W19), 1 86 191. LiDAR A 01.15.01, http://www.airborne1.com/technology/ LiDAR accuracy.shtml ASPR S LiDAR Guidelines (2005), Horizontal accuracy reporting for LiDAR data, http://www.asprs.org/society/committees/standards/Horizontal_Accuracy_R eporti ng_for_Lidar_Data.pdf ASPRS LiDAR Guidelines (2004), Vertical accuracy reporting for LiDAR data, http://www.asprs.org/society/committ ees/lidar/downloads/vertical_accuracy_repo rting_for_lidar_data.pdf Baltsavias, E.P. (1999). Airborne laser scanning: basic relations and formulas. ISPRS Journal of Photogrammetry and Remote Sensing 54 (2 3), 199 214. Besl, P., McKay, N. (1992). A method for registration of 3 D Shap es. Transactions of PAMI 14 (2). Besl, P., lain, R. (1986). I nvariant surface characteristics for 3D o bject recognition in r ange i mages. Computer Vision, Graphics, and Image Processing 33(1), 33 80. Caceres J. (2008). Classific ation of building infrastructure and automatic building footprint delineation using Airborne Laser Swath Mapping data PhD Dissertation, University of Florida, Gainesville, USA. Carter, W.E., Shrestha, R.L., Slatton, K.C. ( 2004 ) Photon counting airborne l aser swath mapping (PC ALSM). Proceedings of SPIE, 4th International Asia Pacific Environmental Remote Sensing Symposium 5661, 78 85. Carter, W.E., Shrestha, R.L., Slatton, K.C. ( 2007 ) Geodetic laser scanning. Physics Today 41 46. Chen,Y., Medioni, G. ( 1991). Object modeling by registration of multiple range image. Proceedings of the 1992 IEEE International Conference on Robotics and Automation, 2724 2729. Congalton, R. G. Green, K. (1999). Assessing the accuracy of remotely sensed d ata: principles and practices. Lewis Publications, Boca Raton, FL. PAGE 120 120 Colomina, I. (2004). Modern sensor orientatio n technologies and procedures. International Archives of Photogrammetry and Remote Sensing, ISPRS Conference. Cramer, M. and Stallmann D. (2002). System calibrat ion for direct g eoreferencing. International Archives of Photogrammetry, Remote Sensing and Spat ial Information Sciences, (34) (3A), 79 84. Fernadez, J.C., Singhania, A., Caceres, J,. Slatton, K.C., Starek, M.J., Kumar, R. ( 2007 ) An overview of LiDAR poin t cloud processing software. University of Florida GEM Center: Report_2007_12_001. (available at http://www.ncalm.cive.uh.edu/ as of November 2010 ). Filin, S. (2003). Recovery of systematic biases in laser altime t ry data using natural surfaces. ISPRS Journal of Photogrammetric Engine ering and Remote Sensing, 69(11), 1235 1242. Friess, P., 2006. Toward a rigorous methodology for airborne laser mapping. Proceedings International Calibration and Orientation Workshop EuroCOW2006, Castelldefels, Spain, 25 27 January, 7p. Gelfan d N., Rusinkiewicz, S. (2003). Geometrically stable sampling for the ICP Algorithm. Proc. International Conference on 3 D Digital Imaging and Modeling. 6 10 October, Banff, Canada Ghatzke, T. D ., Grimm, C. M. (2006). Estimating curvature of triangular meshes. International Journal of Shape Modeling 12(1), 1 28. Gruen, A, Akca, D. (2005). Least squares 3D surface and curve matching. ISPRS Journal of Photogrammetry & Remote Sensing 59, 151 174 Ha bib, A. et al. (2006). Automatic Surface Matching for the Registration of LiDAR Data and MR Imagery, ETRI Journal 28 ( 2 ), 162 174 Habib, A. et al. (2005). Photogrammetric and LiDAR data regis t ration using linear features. Photogrammetric Engineering & Rem ote Sen s ing 71 (6), 699 707 Huisi ng, E. J., Pereira, L. M. G. ( 1998 ) Errors and a ccuracy e stimates of l aser d ata a cquired by various laser s canning systems for t opographic a pplications, ISPRS J. of Photogrammetry and Remote Sensing, 53(5), 245 261. Kilian J. Haala, N., Englich M. (1996). Capture and evaluation of airborne laser scanner data. International 208 The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences. Vol. XXXVII. Part B1. Beijing 2008 Archives of Photogrammetry and Remote Sensing 31(B3), 383 388. PAGE 121 121 Latypov, D. (2002). Estimating r elative LiDAR a ccuracy i nformation from o verlapping f light lines. ISPRS Journal of Photogrammetric Enginee ring and Remote Sensing 56(4), 236 245. Lee, J., Yu, K., Kim, Y. an d Habib, A. (2007) Adjustment of discrepancies b etween LiDAR d ata s trips using l inear features Luzum, B. J., Slatton, K. C., Shrestha, R. L. (2005). Analysis of s patial and temporal s tability of a irborne l aser s wath m apping d ata in f eature s pace", IEEE Transactions in Geoscience and Remote Sensing 43 ( 6 ) 1403 1420. Maas, H. G. (2000). Least s quares m atching with a irborne l aser s canning d ata in a TIN s tructure International Archives of Pho togrammetry and Remote Sensing 33(B3/1), 548 555. Pfeifer, N., Elb erink, S. O., Filin S. (2005). Automatic t ie e lements d etection for l aser s canner s trip a djustment. International Archives of Photogrammetry and Remote Sensing, 36(3/W3), 1 682 1750. Rusi nkiewicz, S., Levoy, M. (2001). Efficien t variants of the ICP Algorit hm. Proc. of 3D Digital Imaging and Modeling ( 3DIM 2001 ). 28 May 1 June Quebec City, Canada. Sampath, A., Shan, J. ( 2007 ) Building boundary tracing and regularization from airborne LiDAR point clouds. P hotogrammetric Engineering and Remote Sensing 73 (7), 805 812. Schenk, T. (2001). Modeling and analyzing systematic err ors in airborne laser scanners. Technical Notes in Photogrammetry (19). The Ohio State University, Columbus, USA. Shrestha, R.L., Carter, W.E., Lee, M., Finer, P., and Sartori, M. ( 1999 ) Airborne l aser s wath m apping: accuracy assessment for surv eying and mapping applications. Journal of American Co ngress on Surveying and Mapping 59 (2), 83 94. Shrestha, R.L., Carter, W.E., Sartori, M., Luzum, B.J., and Slatton, K.C. (2005). Airborne la ser swath mapping: quantifying changes in sandy beaches over time scales of week s to years. ISPRS Journal of Photogrammetry and Remote Sensin g 59 (4), 222 232. Skaloud, J. and D. Lichti (2006). Rigorous approach to bore sight self calibration in airborne laser scanning. ISPRS Journal of Phot ogrammetry and Remote Sensing 61(1), 47 59. Starek, M., Vemula, R., Slatton, K.C. Shrestha, R.L., Carter, W.E. (2007). Automatic f eature e xtraction from a irborne LiDAR m easurements to i dentify c ross s hore m orphologies i ndicative of beach erosion. IEEE Geoscience and Remote Sensing Symp. (IGARSS) PAGE 122 122 Slatton, K., C., Carter, W. E., Shrestha, R. L., Dietrich, W. (2007). Airborne Laser Swath Mapping: Achieving the resolution and accuracy required for geosurficial G eophys ical Research Letters 34, L23S10, doi : 10.1029 /2007GL031939 Triglav Cekada, M., Crosilla, F., and Kosmatin Fras, M. (2009). A simplified analytical model for a priori LiDAR point positioning error estimation and a review of LiDAR error sources. Pho togrammetric Engineering and Remote Sensing 75( 12 ) 1425 1439 Wehr, A., Lohr, U. (1999). Airborne Laser Scanning an Introduction and Overview, ISPRS Journal of Photogrammetry & Remote Sensing 54 ( 2 3 ), 68 82. PAGE 123 123 BIOGRAPHICAL SKETCH Pravesh Kumari was born and raised in a northern state of India called Himachal Pradesh and is the eldest of two children She graduated from Govt. Senior Secondary School, Kangoo. She attended the National Institute of Technology, Hamirpur, where s he graduated with a B. Tech. Honors in Civil Engineering in 2002. She received her M. Tech. in Geoinformatics from Indian Institute of Technology (IIT) Kanpur India in 2004. She completed h er research work for her M. Tech. degree at the University of Stuttgart, Germany through t he student exchange program from 2003 2004 In 2004, Pravesh applied to the Geosensing S ystems E ngineering Ph.D. program at the University of Florida to pursue h er interest in remote sensing. Upon completion of h er Ph.D., Pravesh plans to pursue a career i n academia. Pravesh painting badminton and reading books 