Copyright by Pierre Gueudet 2004 The Influence of Post-Spacing Density of DEMs Derived from LIDAR on Flood Modeling by Pierre Gueudet, B.S Thesis Presented to the Faculty of the Graduate School of The University of Texas at Austin in Partial Fulfillment of the Requirements for the Degree of Master of Science in Engineering The University of Texas at Austin May 2004 The Influence of Post-Spacing Density of DEMs Derived from LIDAR on Flood Modeling Approved by Supervising Committee: Dedication I dedicate this work to my family and to B?n?dicte for their unconditional support and to my African friends from the village of Tansilla in Burkina-Faso who helped me to give a sense to my work as an engineer and to bring my stone in the construction of this world. iv Acknowledgements A portion of this research was funded under a NASA grant (NAG13-02049). I want to acknowledge Jean Parcher (USGS), Bruce Davis (NASA), Judith Berglund (Lockheed), Lee Estep (Lockheed) for their participation in this project. Special thanks to Dr. Gordon Wells and my advisor Dr. David Maidment for their assistance, encouragement, and for their great help throughout the research work. I would also like to thank Amy Neuenschwander for her technical support and guidance as well as Alicia Fogg who spent a lot of time correcting the first draft of my thesis and who gave me very useful feedback. Finally, thanks to my friends at Pickle Research Center for making CRWR a fun place to work. March 30, 2004 v Abstract The Influence of Post-Spacing Density of DEMs Derived from LIDAR on Flood Modeling Pierre Gueudet, M.S.E The University of Texas at Austin, 2004 Supervisor: David R. Maidment The primary objective of the research is to determine the optimal post-spacing for LIDAR-derived digital elevation models (DEMs) that is required to achieve different levels of accuracy in the prediction of flood risk using hydraulic models. For the study, high spatial resolution LIDAR data were collected by the University of Texas Airborne Laser Terrain Mapping System and decimated to generate a variety of DEMs with different resolutions for test and evaluation. The data were entered as input to FEMA-approved hydraulic models at varying resolutions to determine the sensitivity of the models to the changes in the densities of the LIDAR ground elevation points. Data were collected in Brownsville, Texas, in the area of the North Main Drain. A flood model was vi developed using HEC-RAS to delineate the 2-year, 5-year and 25-year floodplains of this drain. For each of the floodplains, 70 simulations were run using different densities of LIDAR-derived ground elevation points. These varying datasets density were also used to delineate watersheds. The work flow sequence needed to produce the flood maps was automated using the ESRI ArcGIS 9 Model Builder. Concerning floodplain delineation, results reveal that below a certain density (5 ground points per 100m 2 which corresponds to a 4.34 meter post-spacing), the LIDAR data become problematical for use in the creation of accurate flood hazard maps. At lower ground point densities, flow obstructions appear on the 3D cross-sections derived from LIDAR, and the inundation polygons expand to unrealistic proportions. With densities greater than 5 ground points per 100m 2 , the influence of increasing LIDAR point density is weak, and the accuracy of the floodplain resulting from a simulation depends upon the initial conditions of the model run. Concerning watershed delineation, no trend was observed as the density of LIDAR points decreases. This result could be explained by the fact that the area of study is extremely flat. vii Table of Contents List of Tables ...........................................................................................................x List of Figures........................................................................................................ xi CHAPTER ONE : INTRODUCTION.....................................................................1 1.1 Background...............................................................................................1 1.2 Motivation and objectives.........................................................................2 1.3 Report structure.........................................................................................3 CHAPTER TWO : STATE OF THE ART OF LIDAR TECHNOLOGY ..............4 2.1 Overview of lidar technology ...................................................................4 2.1.1 Principle of LIDAR.......................................................................5 2.1.2 Description of a LIDAR system ...................................................7 2.1.3 Data classification.........................................................................8 2.2 Review of LIDAR applications ..............................................................12 2.2.1 The main advantages of LIDAR.................................................12 2.2.2 The different applications of LIDAR..........................................13 CHAPTER THREE : METHODOLOGY .............................................................16 3.1 Overview of the methodology ................................................................16 3.1.1 Data collection equipment ..........................................................18 3.1.2 LIDAR data classification...........................................................18 3.1.3 Decimation..................................................................................19 3.1.4 Interpolation to a DEM ...............................................................23 3.1.5 Extraction of the 3D Cross sections and delineation of the watersheds...................................................................................30 3.1.6 Simulation and visualization.......................................................32 3.1.7 Analysis of the results in ArcGIS ...............................................34 3.2 Automation with the ESRI ArcGIS Model Builder................................35 3.2.1 Creation of a tool to decimate LIDAR data and a tool to import geometric data into HEC RAS....................................................37 3.2.2 Description of the model.............................................................39 viii 3.2.2 Looping with the Model Builder.................................................44 CHAPTER FOUR : CONSTRUCTION OF THE HYDRAULIC FLOOD MODEL .......................................................................................................................48 4.1 Description of the area of LIDAR survey...............................................48 4.1.1 The hydrologic landscape in Brownsville...................................48 4.1.2 The hydrologic landscape in Matamoros....................................51 4.1.3 Climate and historical floods ......................................................52 4.1.4 The North Main Drain ................................................................53 4.2 One dimensional hydraulic modeling theory..........................................54 4.2.1 Energy equation ..........................................................................55 4.2.2 Boundary conditions ...................................................................57 4.2 HEC-RAS hydraulic model ....................................................................58 4.2.1 HEC-RAS Geometric data..........................................................58 4.2.2 Calibration of the model .............................................................65 CHAPTER FIVE : RESULTS AND DISCUSSION.............................................69 5.1 Comparison of floodplain and watershed delineation as the density of LIDAR points decreases ......................................................................69 5.1.1 Floodplain results........................................................................70 5.1.2 Results concerning watersheds ...................................................76 5.2 Comparison of the number of buildings flooded as the density of LIDAR points decreases ...................................................................................77 5.3 Computing resources ..............................................................................80 CHAPTER SIX : CONCLUSIONS AND RECOMMENDATIONS ...................82 6.1 Conclusions.............................................................................................82 6.2 Recommendations...................................................................................83 APPENDIX A: REULTS.......................................................................................85 A.1 Watersheds delineation results (to save space, only a portion of the results are presented).......................................................................................85 A.2 Two-year flood results (to save space, only a portion of the results are presented).............................................................................................87 A.3 Five-year flood results (to save space, only a portion of the results are presented).............................................................................................89 ix A.4 Twenty five-year flood results (to save space, only a portion of the results are presented).......................................................................................91 APPENDIX B: CODE ...........................................................................................92 B.1 Decimate tool code.................................................................................92 B.2 Run HEC-RAS tool code .......................................................................93 B.3 Looping Tool code .................................................................................96 GLOSSARY ........................................................................................................102 BIBLIOGRAPHY................................................................................................103 VITA....................................................................................................................105 x List of Tables Table 3-1: Value of the cell size of the DEMs used to interpolate the LIDAR ground points......................................................................................................................... 46 Table 3-2: Boundary conditions for the 2-year 5-year and 25-year flood ........................ 47 Table 3-2: Boundary conditions (flow and water surface elevation at the downstream cross section of the North Main Drain) for the 2-year, 5-year and 25-year flood .... 67 xi List of Figures Figure 1-1: Impact of the resolution of topographical data on the accuracy of the delineation of a floodplain using a hydraulic model. The same simulation was run with high resolution data (left picture) and low resolution topographical data (right picture). .......................................................................................................................2 Figure 2-1: Illustration of the principles of data collection with LIDAR (Maidment, 2001) ..................................................................................................................................... 6 Figure 2-2: A typical airborne laser scanning system (A. Wehr and U.Lohr 1999)........... 8 Figure 2-3: DEMs obtained before LIDAR data classification (upper left), after the vegetation removal (upper right), and after the buildings and vegetation removal (lower left). ............................................................................................................... 11 Figure 3-1: Summary of the methodology used to do the sensitivity analysis ................. 17 Figure 3-2: LIDAR points in text format. From left to right, the column represents the point ID, and the x-y-z coordinate in feet................................................................. 19 Figure 3-3: Production of 70 datasets of LIDAR ground points at varying densities from 104 points per 100m 2 to 1.5 points per 100m 2 .......................................................... 20 Figure 3-4.a: Influence of the velocity of the aircraft on the density of LIDAR points ... 20 Figure 3-4.b: Influence of the altitude of the aircraft on the density of LIDAR points.... 21 Figure 3-5: Illustration of directional spatial decimation, assuming perfectly aligned LIDAR points............................................................................................................ 22 Figure 3-6: Different densities of LIDAR ground points obtained................................... 23 Figure 3-7: Interpolation of LIDAR ground points in 70 DEMs of different resolution.. 24 Figure 3-8: Interpolation under bridges without barrier. After interpolation, bridges are still present even if LIDAR points corresponding to bridges were removed before. 27 Figure 3-9: Interpolation under bridges with barriers....................................................... 29 Figure 3-10: Extraction of the 3D cross sections of the North Main Drain from the DEM by intersecting the 2D line with the edges of the cells of the DEM ......................... 31 Figure 3-11: Delineation of watersheds with different resolution DEMs. The boundaries of the DEMs are symbolized by the square. ............................................................. 32 Figure 3-12: Flood polygon resulting of the HEC RAS simulation ................................. 33 Figure 3-13: Visualization of the 25-year floodplain in ArcGIS, with location of buildings.................................................................................................................... 34 Figure 3-14: Illustration of the definition of the error ...................................................... 35 Figure 3-15: Example of a flow diagram created with the Model Builder....................... 36 Figure 3-16: Interface of the tool "Decimate" .................................................................. 38 Figure 3-17: Interface of the tool "RunHecRas"............................................................... 39 Figure 3-18: Flow diagram of the automated actions required to decimate and interpolate the LIDAR points...................................................................................................... 40 Figure 3-19: Extraction of the 3D cross sections of the North Main Drain from a DEM 41 Figure 3-20: Creation of the floodplain by intersecting the flood polygons produced by HEC RAS and the DEMs.......................................................................................... 42 Figure 3-21: Delineation of the watersheds of the area of study with the DEMs corresponding to different densities of LIDAR ground points. ................................ 43 xii Figure 3-22: Model used to run the sensitivity analysis. The inputs are the LIDAR ground points. The outputs are: the watershed, the buildings flooded and the floodplains. The key parameters (P) are: the boundary conditions, the number n and the cell size. ................................................................................................................................... 45 Figure 4-1: Satellite image (SPOT) of Brownsville and Matamoros ............................... 49 Figure 4-2: Main resacas and drains in Brownsville with their watersheds ..................... 50 Figure 4-3: North Main Drain near Boca Chica Boulevard.............................................. 54 Figure 4-4: Cross section of the North Main Drain subdivided into units with different n values ........................................................................................................................ 57 Figure 4-5: The river system schematic (stream banks, stream flowpath and stream centerline) and their table of contents....................................................................... 60 Figure 4-6: A total of 48 cross sections were used to compute the flow of the North Main Drain. ........................................................................................................................ 61 Figure 4-7: Shapefile of landuse with the cross sections of the North Main Drain.......... 63 Figure 4-8: Table of contents of the laduse shapefile used to extract the Manning's n values ........................................................................................................................ 64 Figure 4-9: Calibration of the hydraulic flood model by adjusting the flow at the downstream cross section so that the computed water surface profile fits with the displayed water surface profile. ................................................................................ 67 Figure 5-1: Illustration of the definition of the error ........................................................ 69 Figure 5-2: Example of the differences in location and extent between the 25-year floodplain computed with highest (left) and lowest (right) density of ground LIDAR points......................................................................................................................... 70 Figure 5-3: Calculation of the error in the delineation of the flood polygons as the density of LIDAR points changes for the 2-year flood ......................................................... 71 Figure 5-4: Calculation of the error in the delineation of the flood polygons as the density of LIDAR points changes for the 5-year flood ......................................................... 72 Figure 5-5: Calculation of the error in the delineation of the flood polygons as the density of LIDAR points changes for the 25-year flood ....................................................... 73 Figure 5-6: Apparition of flow obstructions for low densities of LIDAR points ............. 74 Figure 5-7: The appearance of flow obstructions on the 3D cross sections when using a low density of LIDAR ground points creates discontinuities in the computed water surface profiles due to artificial flow obstruction..................................................... 75 Figure 5-8: Calculation of the error in the delineation of watersheds as the density of LIDAR points changes.............................................................................................. 76 Figure 5-9: Watershed of the North Main Drain computed with the highest density of LIDAR ground points (104 points per 100m 2 , left) and the lowest density of LIDAR ground points (1.5 points per 100m 2 , right).............................................................. 77 Figure 5-10: Calculation of the number of buildings mistakenly flooded by the 2-year flood as the density of LIDAR ground points changes............................................. 78 Figure 5-11: Calculation of the number of buildings mistakenly flooded by the 5-year flood as the density of LIDAR ground points changes............................................. 79 Figure 5-12: Calculation of the number of buildings mistakenly flooded by the 25-year flood as the density of LIDAR ground points changes............................................. 80 xiii Figure 5-13: Evolution of the computation time required to decimate the LIDAR ground points and compute the floodplain as density of LIDAR ground points decreases.. 81 1 CHAPTER ONE: INTRODUCTION 1.1 BACKGROUND Many ancient cities were founded along or near the mouths of rivers, since these strategic trade positions brought economic prosperity. However, rivers can also threaten the local population?s security when floods occur. Therefore, early in history, dams, levees and other structures were built in order to fight these natural disasters. In the recent decades, floods have hit homes and businesses all over the world. The effort to control Mother Nature in this battle against floods is now led by engineers. Powerful computers provide efficient tools to describe or predict flooding phenomena. A flood model consists of two parts: the hydrologic model describes the circulation of water from the atmosphere to the river, and the hydraulic model describes the behavior of water within the river channel. The hydraulic model is the primary focus of this thesis. The level of accuracy of a modeled floodplain depends largely on the resolution of the topographical data (Figure 1-1). Therefore, acquiring spatially dense elevation data is critical for large-scale studies. Airborne Light Detection and Ranging (LIDAR) remote sensing has recently become a standard method for the acquisition of high-resolution topographic data. These data are a set of elevation points, spread over the area of study. 2 Figure 1-1: Impact of the resolution of topographical data on the accuracy of the delineation of a floodplain using a hydraulic model. The same simulation was run with high resolution data (left picture) and low resolution topographical data (right picture). 1.2 MOTIVATION AND OBJECTIVES One of the federal government agencies interested in using LIDAR data is the Federal Emergency Management Agency (FEMA), which certifies and maintains the national Flood Insurance Rate Maps (FIRMs). These maps delineate Special Flood Hazard Areas (SFHA), which are areas of land that would be inundated by a flood that has a 1% chance of occurrence in any given year. In 1997, in some cases, FEMA identified a need to update the FIRMs that are currently more than ten years old. They plan to use last return LIDAR data to update their floodplains and need to know the optimum density of points to use. The density of LIDAR data is especially important because it is the greatest cost factor in LIDAR data acquisition / processing, and because increasing the density of points can improve accuracy. A higher density requires one or more of the following factors: a sensor system with a high pulse rate (e.g. 50,000 pulses per second), lower elevation flights (and therefore more flight lines), a lower ground speed, a narrower scan angle, or a combination of these. These techniques significantly increase the cost of 3 LIDAR data acquisition for a given area. Beyond acquisition costs, significantly more computing resources (processor speed, RAM, storage space, etc.) and staff time are required to process higher posting densities during the Digital Elevation Model (DEM) creation process. How does one determine a balance between LIDAR data product accuracy and project cost? The objective of this study is to conduct a sensitivity analysis to define the optimum density of topographical data and thus, the best expenditure strategy required to achieve different levels of accuracy in the prediction of flood risk using hydraulic models. Currently there are few guidelines to answer this question. Some studies (North Carolina, 2002; Hodgson, et al., 2003) have attempted to document the accuracy of Digital Terrain Models (DTMs) derived from LIDAR. Fewer studies have examined the effect of DTM accuracy on hydraulic modeling (Kenward, et al., 2002). No study on the effects of LIDAR data density on flood mapping was found. 1.3 REPORT STRUCTURE Following this introduction, Chapter 2 provides more detailed information about LIDAR technology to justify the methodology, which is presented in Chapter 3. Chapter 4 is devoted to the construction of a flood model using the HEC River Analysis System (HEC RAS) hydraulic model, which was developed by the US Army Corps of Engineers. The results of the sensitivity analysis are presented in Chapter 5. Finally, Chapter 6 provides concluding thoughts on the use of LIDAR data for mapping and offers recommendations for future research. 4 CHAPTER TWO: STATE OF THE ART OF LIDAR TECHNOLOGY The idea of using an airborne laser for measuring elevations was first conceived in the late 1960?s. However, the actual developments in the technology of airborne altimetric LIDAR (Light Detection And Ranging) only emerged with some parallel advancements in the field of position and attitude measurement using GPS (Global Positioning System) and IRS (Inertial Reference System). Despite its comparatively recent emergence, the pace of development in LIDAR has been very rapid. The technology is now considered the most appropriate method for altimetry. This chapter briefly describes the principles and technical issues of this technology. Furthermore, the varied successful applications of this technology are reviewed. 2.1 OVERVIEW OF LIDAR TECHNOLOGY Topographic data are the core of many GIS (Geographical Information System) projects. The accuracy and functionality of many GIS projects rely to a large extent on the accuracy of topographic data and the efficiency with which they can be collected. Furthermore, the data collection consumes a major slice of the project resources in terms of both time and money. Topographic data collection, therefore, is a significant part of any GIS project. The conventional methods of topographic data collection include land surveying and aerial photogrammetry. More recently, attempts have been made to use satellite stereogrammetry for this purpose. However, all the above techniques are limited in their accuracy, cost-effectiveness, time-consumption, feasibility, and applicability. The recently emerged technique of airborne altimetric LIDAR has gained considerable acceptance in both scientific and commercial communities as a tool for topographic 5 measurement. This technique has the potential to remove several bottlenecks imposed by earlier methods. 2.1.1 Principle of LIDAR The basic concepts of an airborne LIDAR mapping system are simple. The airborne LIDAR instrument transmits the laser pulses while scanning a swath of terrain, usually centered on and co-linear with the flight path of the aircraft carrying the instrument. The scan direction is orthogonal to the flight path. The round trip travel times of the laser pulses from the aircraft to the ground are measured with a precise interval timer. Since the speed of light is known, the time intervals can be converted into range measurements, i.e., the distance from the LIDAR instrument to the ground point struck by the laser pulse, employing the velocity of light. The position of the aircraft at the instant of firing the pulse is determined by differential GPS. Rotational positions of the laser pulse direction are combined with aircraft roll, pitch, and heading values determined with an inertial navigation system (INS) and with the range measurements to obtain range vectors from the aircraft to the ground points. When these vectors are combined with the aircraft locations, they yield accurate coordinates of the points on the surface of the terrain. Figure 2-1 illustrates the principle of LIDAR. 6 Figure 2-1: Illustration of the principles of data collection with LIDAR (Maidment, 2001) A typical LIDAR campaign involves the following steps: ? Flight planning, i.e. setting LIDAR instrument and aerial platform parameters to control the density and coverage of topographic measurements ? Defining ground control points (GCPs) to place reference receivers for differential GPS positioning ? Instrument calibration before, during, and after the flight to ensure accuracy of data collected ? Data collection, i.e. obtaining INS, GPS, laser range and scan measurements GPS ground reference station Aircraft GPS GPS satellites Flight direction L a s e r s c a n l i n e s ALTM laser and IMU 7 ? Data processing to determine aerial platform location using GPS and INS measurements and using the laser range and scan measurement to yield triplets (i.e., x, y, and z) for each ground point struck by the laser in World Geodetic System (WGS-84) ? Quality assurance to determine and quantify the errors present in data and, if needed, elimination or minimization of the errors ? Generation of data products, i.e., DTM, DEM, contour plots and 3D visualization 2.1.2 Description of a LIDAR system All LIDAR systems use some means to measure the distance between the sensor and the illuminated spot on ground. As shown in Figure 2-2, a typical laser scanner can be subdivided into the following key units: laser ranging unit, opto-mechanical scanner, controlling and processing unit, and position and orientation system (POS). The ranging unit comprises the emitting laser and the electro-optical receiver. The transmitting and receiving apertures, typically 8?15 cm in diameter, are mounted so that the transmitting and receiving paths share the same optical path. This assures that object surface points illuminated by the laser are always in the field of view (FOV) of the optical receiver. The narrow divergence of the laser beam defines the instantaneous field of view (IFOV). Typically, the IFOV ranges from 0.3 milliradians to 2 milliradians. The theoretical physical limit of the IFOV is determined by diffraction of light, which causes image blurring. The position and orientation system POS is obtained by an integrated differential GPS (DGPS) and an inertial measurement unit IMU. 8 Figure 2-2: A typical airborne laser scanning system (A. Wehr and U.Lohr 1999). 2.1.3 Data classification Based on an integrated multi-sensor system with GPS and INS, an airborne laser terrain mapper provides the 3D position of each laser beam spot on the Earth's surface. The most important feature of recent airborne LIDAR systems is their ability to differentiate the first pulse reflections from the last pulse reflections. A laser pulse that is fired over vegetation usually has multiple reflections. Some components of the laser pulse may be reflected by leaves or branches of trees often represented in the first returning pulse. Others may be reflected by the ground, and the last returning pulse is most likely to be reflected by the terrain surface beneath trees. This ability to differentiate first and last pulse reflections, coupled with specific interpolation and data filtration (signal processing) methods, allows one to separate the LIDAR ground points from the points on buildings and vegetation. Then, a DTM can be computed. However, no general 9 algorithm exists to extract the buildings or the vegetation: different algorithms are used to classify the data for each feature class. The principle, however, is always used and is described below. For classification of buildings, first and last pulse data are expected to be similar in return time and intensity because of the mostly solid surface of buildings. Along 3D discontinuities, it is known that some first pulses represent reflections from the roof surface, while the last pulse may result from a reflection on the ground. Another characteristic of roofs in laser images is that they show up with a more or less constant surface orientation, which means that gradients of the range data should be homogenous in the local area covered by the roof. Finally, the typical size of buildings is a good property to distinguish the buildings from other mostly smaller 3D objects. Therefore, based on the characteristics described above, interpolation methods can be used to extract the points corresponding to buildings. As already mentioned, a key feature for the detection of trees is the ability of the laser to penetrate vegetation to a certain extent. First-return pulses are more likely to reflect off leaves or branches of trees, while last return pulses often reflect off the ground below the trees or forest. Similar to the well-known Normalized Difference Vegetation Index (NDVI) in multispectral image data classification, which employs near infrared and red image channels, a corresponding vegetation index can be defined based on the first and last pulse images (Equation 1). ND is the Normalized Difference of first (FR) and last (LR) return travel time. A high ND indicates the presence of vegetation. The Normalized Difference image is an LRFR LRFR ND + ? = (1) 10 index indicating vegetation. This equation is the basis for vegetation detection and extraction. Figure 2-3 gives an example of a DEM produced before data classification, after removal of the buildings, and after removal of the buildings and the vegetation. 11 Figure 2-3: DEMs obtained before LIDAR data classification (upper left), after the vegetation removal (upper right), and after the buildings and vegetation removal (lower left). 12 2.2 REVIEW OF LIDAR APPLICATIONS This section describes the main advantages of LIDAR technology compared with traditional topographical data collection methods. A short overview of the applications of LIDAR is given. 2.2.1 The main advantages of LIDAR These are the main advantages of LIDAR technology compared with the traditional methods of topographic data. ? Collection accuracy: accuracy on the order of 10 to 15 cm in the vertical direction and 50 to 100 cm in the horizontal direction is claimed by manufacturers and has been demonstrated by many field studies. ? Time of data acquisition and processing: The data capture and processing time is significantly less for LIDAR compared to other techniques. LIDAR can allow surveying rates of up to 90 km 2 per hour (Environment Agency, 1997), with post-processing times of two to three hours for every hour of recorded flight data (Martin and Gutelius, 1997). ? Additional data: Besides relief information, laser reflectance may be used to generate intensity images to help in classifying the terrain features (Schreier et al., 1985). Further, the systems can have fluorosensors, allowing pollutant identification and chlorophyll mapping (Environment Agency, 1997; Measures, 1984). ? Canopy penetration: Unlike photogrammetry, LIDAR can see below the canopy in forested areas and provide topographic measurements of the surface underneath. Additionally, LIDAR generates multiple returns from single pulse travel, thus providing information about the vegetation understory. 13 ? Data density: LIDAR can measure subtle changes in terrain, as it generates a very high data density (due to firing of 2000 - 80000 pulses per second). ? Ground control point independence: Each LIDAR pulse is individually georeferenced using the onboard GPS, INS, and laser measurements. Only one or two GPS ground stations are required for improving the GPS accuracy by the differential method. Independence from GCPs makes LIDAR an ideal method for inaccessible or featureless areas like wastelands, ice sheets, deserts, forests, and tidal flats. ? Digital compatibility: Data produced from LIDAR flights are in digital format with Easting, Northing, and Altitude values of each laser target, which makes importing of data to GIS and other image processing packages straightforward. ? Cost: One of the major hindrances in the use of LIDAR had been the cost of the equipment. However, in recent years the purchase price of these instruments has been reduced so that cost is no longer a barrier to companies capable of investing in standard aerial photogrammetry equipment (Martin and Gutelius, 1997). Furthermore, with more and more users opting for LIDAR, the cost of the system and operation is likely to go further down. On the basis of overall performance evaluation of available topographic techniques for coastal terrain Mason et al. (1999), found that LIDAR could achieve good performance at a lower cost. 2.2.2 The different applications of LIDAR A literature review reveals that LIDAR has been used in gathering information to manage various natural and artificial disasters. In addition, due to its unique data collection properties, LIDAR is also finding use in several other application areas not 14 deemed feasible until now. The following paragraphs describe the state-of-the-art in LIDAR technology, its advantages and disadvantages, and an overview of its present and prospective uses in disaster management. ? Floods: High-resolution and accurate LIDAR data are suitable for improving the performance of flood models by providing a more reliable initial boundary condition (Bates, 1999). LIDAR data having multiple returns help in generating and understanding the 3D structure of obstructions (i.e. surface roughness, vegetation, buildings, and other structures). This information can yield the friction coefficient over the various parts of a floodplain (Cobby, 1999). LIDAR data are also being employed for flood hazard zoning (Hill et al., 2000). FEMA is using LIDAR on a mandatory basis to create Digital Flood Insurance Rating Maps. ? Coastal applications: LIDAR has generated considerable interest among coastal researchers as a topographic tool. Highly accurate, dense, and rapidly obtained data sets are most suitable for coastal applications like sediment transport, coastal erosion, and coastal flood models Brock et al., 1997; Gutierrez et al., 1998). ? Bathymetry: LIDAR in its bathymetric form (shorter wave length) can map the bed topography up to a depth of 70 m in extremely clear water (Wehr and Lohr, 1999). This information is useful for determining the siltation on navigation canals and ports and for planning the construction details. ? Hydrology: LIDAR data can be used to quantify gully and stream channel cross sections and roughness, gully and stream bank erosion and channel 15 degradation, to estimate soil loss from gully or channel banks and to measure channel and floodplain roughness and cross sections for estimating flow rates (Ritchie, 1996). Further, coastal channels, which are difficult to map otherwise, can be automatically quantified using LIDAR data (Lohani, 1999). ? Glaciers and avalanches: LIDAR has been found ideally suited for mapping glacial topography (Krabill et al., 1995; Kennett and Eiken, 1997) and ice sheet velocities (Abdalati and Krabill, 1999). The above studies have amply shown that LIDAR can be used to monitor glacier movement and snow accumulation and predict the onset of avalanches. The data set can further be employed to estimate the risk from a particular avalanche (Wehr and Lohr, 1999). ? Landslides: LIDAR has made it possible to monitor and predict slope failure by rapidly obtaining highly accurate and dense elevation data. In post-slide conditions, rapid damage assessment and mapping can be realized using LIDAR. ? Forest mapping: LIDAR?s unique feature of producing multiple returns from the canopy top, understory, and the ground has attracted many to use it for estimating forest biomass, timber volume, and other parameters (Nelson et al., 1984). ? Volcano monitoring: Ridgway et al. (1997) have shown that subtle systematic changes (uplift of up to 4 cm per year) in volcano dome height can be monitored from time to time using LIDAR. 16 CHAPTER THREE: METHODOLOGY 3.1 OVERVIEW OF THE METHODOLOGY The objective of this project is to define the optimum density of LIDAR ground points, and thus the optimum cost required to achieve different levels of accuracy in the prediction of flood risk using hydraulic models. Ideally, a separate LIDAR dataset should be acquired for each density of ground points under investigation, which would involve holding as many sensor variables constant as possible (such as scan angle and pulse rate) and changing only one variable (such as altitude or velocity) as needed. However, this project takes advantage of an existing high density LIDAR dataset for the study area (approximately one ground point per square meter). The original high resolution LIDAR dataset was decimated to generate 70 datasets of ground points with a resolution varying from 104 points per 100m 2 to 1.5 points per 100m 2 . These datasets were used to generate a variety of DEMs. The cell size of the DEM was defined so that there is an average of one point per cell. The DEMs of varying resolutions were entered as input to a FEMA-approved hydraulic model to determine the model?s sensitivity. The extent and the location of the floodplains obtained with the different LIDAR ground point datasets as well as the number of buildings flooded were compared. One ?side product? of the DEMs is the watershed delineation that can easily be derived from a DEM. Therefore, even if watersheds are required only for hydrologic modeling, the DEMs were used to examine the influence of the LIDAR point density on watershed delineation by comparing their extent and their location. Research for this project was conducted in Brownsville, Texas, on an area of the North Main Drain. The following flow diagram 17 (Figure 3-1) summarizes the methodology followed to make the sensitivity analysis. Each step is explained in more detail in the next paragraphs. Figure 3-1: Summary of the methodology used to do the sensitivity analysis 18 3.1.1 Data collection equipment To collect the LIDAR data over the region of metropolitan Brownsville, Texas, an Airborne Laser Terrain Mapping (ALTM) 1225 system was used. The instrument is a scanning LIDAR system developed by Optech, Inc., of Toronto. The ALTM 1225 has the following specifications: ? Operating altitude 410-2,000 m above ground level (AGL) ? Range resolution: 1 cm ? Laser pulse repetition rate: 25 kHz ? Laser scan angle: variable from 0 to 20 degrees from nadir ? Laser scanning frequency variable, 28Hz at the 20 degree scan angle ? Beam divergence variable, 0.2 or 1.0 milliradian ? Simultaneous recording of first return, last return, and intensity of laser reflections 3.1.2 LIDAR data classification The first step in the process is the data classification. The LIDAR points corresponding to vegetation and structures were removed by the Center for Space Research in Austin, TX, to retain only the LIDAR ground points. The LIDAR data are stored as a text file with five columns that contain the x, y and z point coordinates and the intensities of the first and last return of the laser pulses (Figure 3-2). For this project, the geographic coordinate system of the LIDAR data is the North American Datum of 1983, and the projected coordinate system is NAD 1983 State Plane Texas South, FIPS 4205, Feet. 19 Figure 3-2: LIDAR points in text format. From left to right, the column represents the point ID, and the x-y-z coordinate in feet. 3.1.3 Decimation The objective of the decimation process is to use the original high density LIDAR ground points to generate datasets of LIDAR ground points with different densities. Seventy datasets with different densities varying from 104 points per 100m 2 to 1.5 points per 100m 2 were produced (Figure 3-3). 20 Figure 3-3: Production of 70 datasets of LIDAR ground points at varying densities from 104 points per 100m 2 to 1.5 points per 100m 2 A higher density of ground points requires a lower flight elevation (and therefore more flight lines), a lower ground velocity, a narrower scan angle, or a combination of these variables. Changing the altitude or the velocity at which points are sampled has different impacts on the spatial distribution of points. Increasing the velocity of the aircraft affects the density of LIDAR ground points by increasing the distance between the laser scan lines (multiplying the velocity of the aircraft by two multiplies the distance between scan lines by two). Alternately, increasing the altitude of the aircraft increases the distance between two consecutives LIDAR ground points (multiplying the altitude by two multiplies the distance between two consecutive points by two). Figures 3-4.a and 3- 4.b illustrates the concept. Figure 3-4.a: Influence of the velocity of the aircraft on the density of LIDAR points Speed 1 Speed 2 Speed 1 < Speed 2 21 Figure 3-4.b: Influence of the altitude of the aircraft on the density of LIDAR points Given the huge size of the original LIDAR data (10,213,012 LIDAR ground points for the study area, which represents several Gigabits), the decimation process was handled in Microsoft Access for more efficiency. Two common decimation methods exist: random decimation and directional spatial decimation. A previous study (Chow 1997) shows that ?no best decimation method can be suggested.? However, for this study, after different trials, the directional spatial decimation method was preferred instead of the random decimation. ? Random decimation method: in random sampling, each LIDAR point is assigned a random number. They are then sorted in ascending order and the first n numbers of the sorted features are selected, depending on the desired sample size. This process does not simulate any of the physical cases, such as the effect of flight altitude or velocity described above. Indeed, the data are collected with a certain order and regularity and a random decimation does not conserve this regularity: this decimation process creates some regions with no points whereas other regions still have a high density of points. This phenomenon can be observed Altitude 1 Altitude 2 Altitide 1 < Altitude 2 22 especially when a high percentage of points are deleted. Therefore, this method was eliminated. ? Directional spatial decimation method: the original high resolution LIDAR points are classified in the order of collection by the sensor and then decimated by deleting n points every n+1 points (Figure 3-5). This decimation method simulates how the density of LIDAR points decreases as the altitude of the sensor system increases. Seventy sets of LIDAR ground points were generated (Figure 3-6). Figure 3-5: Illustration of directional spatial decimation, assuming perfectly aligned LIDAR points 23 Figure 3-6: Different densities of LIDAR ground points obtained 3.1.4 Interpolation to a DEM LIDAR points are usually interpolated into a Digital Elevation Model (DEM) for flood modeling. Seventy DEMs were produced from the LIDAR ground points (Figure 3- 7). 0 10 20 30 40 50 60 70 80 90 100 0 10 20 30 40 50 60 70 80 90 100 Percentage of the original LIDAR ground points deleted D ens i t y o f LI D A R groun d po i n t s i n poi n t s per 100 m 2 24 Figure 3-7: Interpolation of LIDAR ground points in 70 DEMs of different resolution The problem is to choose an interpolation method and a cell size corresponding to each density of LIDAR points. The interpolation process has three objectives: ? First, buildings and vegetation have been removed after the classification process; there are no LIDAR points at these locations. Therefore, the data need to be interpolated. ? Second, water present in the rivers during the data collection absorbs the laser pulse due to its wave length: the LIDAR system cannot penetrate water. Therefore, the bottom of the river, which is critical for flood modeling, has to be derived by interpolation. ? Third, hydraulic features such as bridges are an obstacle for the laser; therefore, no LIDAR ground points are collected under bridges. Before interpolating the LIDAR ground points, the points corresponding to the bridges need to be removed. 25 There are two main types of interpolation techniques: deterministic and geostatistical. Deterministic interpolation techniques create surfaces from measured points, based on either the extent of similarity (e.g., Inverse Distance Weighted) or the degree of smoothing (e.g., Radial Basis Functions) of the points. Geostatistical interpolation techniques (e.g., Kriging) utilize the statistical properties of the measured points. Geostatistical techniques quantify the spatial autocorrelation among measured points and account for the spatial configuration of the sample points around the prediction location. Used correctly, geostatistical techniques require experience and involve some subjectivity. In addition, the computation time required to interpolate points is usually very high compared with deterministic interpolation techniques. In general, deterministic interpolation techniques are preferred over geostatistical techniques to interpolate LIDAR points. Based on these factors, deterministic interpolation techniques were used in this study. Deterministic interpolation techniques can be divided into two groups, global and local. Global techniques interpolate points using the entire dataset. Local techniques interpolate points using the measured points within neighborhoods, which are smaller spatial areas within the larger study area: they are more appropriate for the purposes of this study since the objective is to show how details on a DEM affect a flood simulation. Inverse Distance Weighted (IDW) interpolation assumes that things that are close to one another are more alike than those that are farther apart. To predict a value for any unmeasured location, IDW uses the measured values surrounding the prediction location. The measured values closest to the prediction location have more influence on the predicted value than those farther away. Thus, IDW assumes that each measured point has a local influence characterized by its weight that diminishes with distance. 26 Based on the former discussion, the IDW was chosen to interpolate LIDAR ground data. The interpolation was done using ArcGIS and the following parameters were used: ? The weighting of each point is inversely proportional to the distance squared between the predicted location and the measured point. ? The shape of the neighborhood is a disc. The specified shape of the neighborhood restricts how far and where to look for the measured values to be used in the prediction. The searching radius is variable and is defined so that 12 sample points are used to interpolate the prediction location. ? The cell size varies with the density of LIDAR ground points. ? The barriers option specifies the location of bridges. The cell size was chosen so that there is on average one LIDAR ground point per cell. This commonly used criterion defines exactly one cell size for a given density. Thus, as the density of LIDAR ground points decreases, the cell size of the corresponding DEM increases. This has a very important impact on the global computation time of the sensitivity analysis since the cell size is the key parameter which defines the computation time. For example, it takes more than 12 hours to create the DEM for the study area with the highest resolution (104 points per 100m 2 ) and only a few minutes to create the DEM for the same area with the lowest resolution (to 1.5 points per 100m 2 ). It was necessary to use barriers for interpolation. Indeed, the presence of bridges was a source of problems since bridges interrupt the surface continuity. The LIDAR points corresponding to the bridges were deleted manually using a mask whereas the LIDAR points located on the banks just next to a bridge were not since they are ground data. This last category of points has a much higher elevation than the bottom of the river 27 and cannot be used to interpolate the elevation of the ground under the bridge. Figure 3-8 describes the problem. Figure 3-8: Interpolation under bridges without barrier. After interpolation, bridges are still present even if LIDAR points corresponding to bridges were removed before. LIDAR data collection Removal of LIDAR points corresponding to bridges Selection of 12 sample points to interpolate an elevation Interpolation into a DEM 28 Therefore, barriers were used to limit the selected set of sample points used to interpolate z values under the bridge to only those samples on the same side of the barrier as the current processing cell. Separation by a barrier is determined by line-of-sight analysis between each pair of points, which means that topographical separation is not required for two points to be excluded from each other's region of influence. The barriers option significantly extends the computation time required to obtain the DEMs but increases the quality of interpolation. Figure 3-9 illustrates how the sample points were chosen to interpolate the ground under the bridge. 29 Figure 3-9: Interpolation under bridges with barriers. LIDAR data collection Removal of LIDAR points corresponding to bridges and addition of barriers Selection of 12 sample points to interpolate an elevation Interpolation into a DEM 30 3.1.5 Extraction of the 3D Cross sections and delineation of the watersheds. A one-dimensional hydraulic flood model was developed in HEC RAS for a small portion of the North Main Drain in Brownsville, Texas. More explanations and theoretical background concerning the way it was constructed are given in the next chapter. Since the objective of the study is to understand the influence of the density of LIDAR ground points on the model, the only parameter that varies in the simulations is the density of LIDAR points. All the other parameters remained constant. The way the model ?sees? the variation of density of LIDAR points is through the cross sections of the river. The location of the cross section is fixed but their shape changes with the resolution. 2D cut lines were edited by hand, and the 3D cross sections of the North Main Drain were extracted from the DEM by locating the points of intersection of the 2D cut lines with the DEM. Different sets of 3D cross section were extracted from the DEMs and used to run a simulation. Figure 3-10 illustrates this concept: 31 Figure 3-10: Extraction of the 3D cross sections of the North Main Drain from the DEM by intersecting the 2D line with the edges of the cells of the DEM With the DEMs, it was also possible to delineate the watersheds for this part of the North Main Drain. These are not required for hydraulic modeling but are used for hydrologic modeling. The accuracy of the delineation of the watersheds is an important parameter in hydrology. Therefore, the influence of the density of LIDAR data on the delineation of the watersheds was studied. For each of the 70 DEMs, the watersheds were delineated as shown in Figure 3-11. . Cross section 32 Figure 3-11: Delineation of watersheds with different resolution DEMs. The boundaries of the DEMs are symbolized by the square. 3.1.6 Simulation and visualization In order to test the sensitivity of the model to the density of LIDAR points under different conditions, the flood model has been developed to delineate the 2-year, 5-year and 25 year floodplain. For each floodplain, 70 simulations were made with the different DEMs corresponding to the different densities of LIDAR points. The result of a simulation in HEC RAS is a flood polygon as shown by Figure 3-12: 33 Figure 3-12: Flood polygon resulting of the HEC RAS simulation This flood polygon was then exported to ArcGIS and intersected with the DEMs to get the corresponding floodplains. Thus, for each flood event (the 2-year, 5-year and 25 year flood), seventy floodplains corresponding to seventy different densities of LIDAR points were generated and analyzed. Figure 3-12 illustrates the process. One of the government agencies interested in the project is the Federal Emergency Management Agency. Beside the extent and the location of the floodplains, they are also interested in the number of buildings flooded in order to estimate the cost of damages. Therefore, a shapefile describing the location of the buildings was intersected by each floodplain to determine the number of buildings flooded (Figure 3-13). 2242.270 2135.371 2082.934 2043.877 1924.482 1860.756 1837.580 1809.467 1765.357 1727.277 1643.130 1550.141 1500.123 1452.494 1408.067 1369.250 1323.618 1255.878 1196.056 1152.308 1105.280 1020.670 980.2624 887.8959 804.3178 743.0336 711.8319 634.6898 599.7890 553.7707 492.1651 448.3111 409.7674 370.5898 308.2700 197.3756 146.0169 101.2411 NorthMainDrain Plan: Plan 04 11/5/2003 34 Figure 3-13: Visualization of the 25-year floodplain in ArcGIS, with location of buildings. 3.1.7 Analysis of the results in ArcGIS Two types of results were analyzed: the influence of the density of LIDAR points on the extent and the location of the computed floodplains and the influence of the density of the LIDAR points on the number of buildings flooded. ? Comparison of the floodplains: for a given flood event, the floodplains obtained with the different densities of the LIDAR data were compared with the floodplain obtained with the highest density of LIDAR points. To compare these polygons and to take into account the differences in extent and location, the error has been defined as follows (Equation 2, Figure 3- 14): Resaca 35 Figure 3-14: Illustration of the definition of the error The error can vary from 0% to 100%. This formula can be understood by analyzing the two following cases: if the flood polygons have the same area but are located at two different places (intersection null), the error is 100% since the area of the intersection of these two flood polygons is zero. The error is null only if the two flood polygons exactly match. ? Comparison of the number of buildings flooded: based on the same principle as the one used to compare the floodplains, the buildings flooded obtained with the highest resolution of LIDAR points were taken as a reference. In order to be more efficient and to avoid errors, all the steps described previously were automated using the ESRI ArcGIS Model Builder. 3.2 AUTOMATION WITH THE ESRI ARCGIS MODEL BUILDER Model Builder is a system in ArcToolbox that allows geoprocessing tasks to be chained together in a single model. An ArcToolbox model consists of a set of processes. 100 )( )(2)()( (%) ? ???+ = PolygonREFArea PolygonREFPolygonAreaPolygonREFAreaPolygonArea Error (2) Polygon REF Polygon ErrorIntersection 36 Each process contains a tool, input data, and output data. The tool performs a spatial operation on the input data and the output from one tool becomes the input for the next tool in the chain. Model Builder uses a visual approach: it uses a diagram or a flow chart to represent a model, and a graphic icon or a node to represent each process component (Figure 3-15). Figure 3-15: Example of a flow diagram created with the Model Builder. The flow diagrams created with the Model Builder are not only a convenient way to construct and modify spatial models but are also an excellent way to document and share one's models with others. Users can rerun their models using different function parameters, thus enabling them to examine how they perform using different sets of values. Therefore, ArcGIS Model Builder is a good tool to support this sensitivity analysis: a model was developed to take as an input the original ground LIDAR data, process them, run a flood simulation and compare the results. The changing parameter is the density of the LIDAR data or the degree of decimation. 37 The entire sensitivity analysis was automated using the ArcGIS 9 Beta 1 Model Builder in order to delineate the 2-year, 5-year and 25-year floodplains and to run 70 simulations for each one. ArcGIS 9 will be release is June 2004. 3.2.1 Creation of a tool to decimate LIDAR data and a tool to import geometric data into HEC RAS The ArcGIS 9 ArcToolbox allows the user to create custom tools in addition to the standard tools provided with ArcGIS. Script tools send geospatial inputs to a custom script that performs a task outside of the functionality provided by standard ArcGIS tools. The scripts may be written in a variety of scripting languages (e.g. Python, VBScript) and are stored on disk. The ArcToolbox script tool contains a link to the script, as well as input/output parameter definition for the script (i.e., what goes into the script, and what comes out of the script). The script tool retrieves geospatial inputs and then sends those inputs to the script where processing takes place. When the script is finished, it sends an output back to the script tool. A TOOL TO DECIMATE LIDAR DATA A script tool was developed to decimate LIDAR data. The input is a Microsoft Access table containing the LIDAR data and the output is a Microsoft Access table containing the decimated data. The tool named ?Decimate? provides the user with two types of decimation method: a random decimation or a directional spatial decimation. To run the tool the user needs to specify the name of the input and output table as well as the decimation method desired and a parameter n. The tool deletes n-1 points every n points (Figure 3-16). 38 Figure 3-16: Interface of the tool "Decimate" A TOOL TO IMPORT GEOMETRIC DATA INTO HEC RAS A script tool was developed to import geometric data into HEC RAS. The input is an .sdf file containing the geometry of the cross section of the river and the banks. The output is an .sdf file containing the results of a HEC RAS flood simulation. The user has to provide the path of the input .sdf file and the output .sdf file as shown by Figure 3-17: 39 Figure 3-17: Interface of the tool "RunHecRas" 3.2.2 Description of the model This section describes the entire model. The input of the model is the original ground LIDAR dataset, and the outputs are the watershed, the floodplain and the buildings flooded computed with a certain density of LIDAR points. Everything else is a parameter or an intermediate output. Running the model with a given set of parameters produces one watershed, one floodplain and the buildings inside this floodplain for a given density of LIDAR points. In this section, all the parameters remain constant. The next section describes how to select certain parameters and to modify them automatically in order to study the sensitivity of the model to the density of LIDAR points. The first part of the model decimates and interpolates the LIDAR points into a DEM. The following flow diagram (Figure 3-18) represents the succession of actions required to achieve this task. 40 Figure 3-18: Flow diagram of the automated actions required to decimate and interpolate the LIDAR points The tool named ?Decimate? takes a table with all the LIDAR points and creates a new table with fewer points by selecting n point every n+1 points. Based on this table, a feature class of the LIDAR points is created in ArcGIS. Then, the tool named ?IDW? interpolates the points into a grid with a given cell size using a shapefile with the barriers required to interpolate under bridges. All the work is done in the North, American datum of 1983 and UTM Zone 14N projected coordinate system. Since the original elevation was in feet with the State Plane Texas South FIPS 4205 projection, this grid was multiplied by the conversion rate to convert from feet to meters and reprojected. Then, the DEM is used to extract the 3D cross sections with the following tool (Figure 3-19): 41 Figure 3-19: Extraction of the 3D cross sections of the North Main Drain from a DEM The geometric .sdf file containing these 3D cross sections and the banks is created by hand using the ?export to HEC RAS? button of the HEC GeoRAS preprocessor. The DLL corresponding to this button was not available so a tool to do this task could not be created. The second part of the model runs the hydraulic simulation in order to compute the floodplain. The flow diagram is presented in figure 3-20: 42 Figure 3-20: Creation of the floodplain by intersecting the flood polygons produced by HEC RAS and the DEMs. The tool called ?run HecRas? imports the .sdf file into HEC RAS and runs a steady state simulation. The results of the simulation (e.g. the water surface elevation at each cross section) are exported to ArcGIS through an .sdf and a .xml file. Then the tool named ?Unclipped WS TIN? creates a surface of the water elevation and the tool named ?LessThan? intersects this surface with the DEM in order to delineate the floodplain. Finally, the floodplain is converted to a feature class. 43 The DEM is also used to delineate the watershed of the area of study in the following flow diagram (Figure 3-21): Figure 3-21: Delineation of the watersheds of the area of study with the DEMs corresponding to different densities of LIDAR ground points. The watershed delineation is handled in several steps: ? Fill Sinks: at times minor discontinuities occur in the DEM which causes pits to form in the terrain. Sinks are filled so that they are leveled with the surrounding terrain. Only tiny sinks are filled, since large sinks, such as lakes, are real sinks that don?t need to be removed from the DEM. ? Flow Direction: since water flows downhill when unobstructed, this command considers every cell in the DEM grid and decides which of the eight cells surrounding the considered cell offers the steepest downslope gradient. 44 ? Flow Accumulation: eventually, as water flows toward lower elevations, streams will begin to form. The flow accumulation function addresses each cell of the DEM and counts how many upstream cells contribute to flow through the given cell. ? Snap Pour points: snaps pour points to the cells of the highest flow accumulation within a specified distance. ? Watershed: finally, using the stream centerline and the flow accumulation grid, the tool delineates the watersheds in the basin. The girded watersheds are converted into polygons. 3.2.2 Looping with the Model Builder Since the objective is to study the sensitivity of the model to the density of LIDAR points, a script tool was created in order to control the parameters of the model and to run the model 70 times with different values for the parameters n and cell size. These 70 simulations were run for three different types of boundary conditions representing the 2-year, 5-year and 25-year flood. Figure 3-22 shows the entire model with its key parameters for the sensitivity analysis and the main outputs. 45 Figure 3-22: Model used to run the sensitivity analysis. The inputs are the LIDAR ground points. The outputs are: the watershed, the buildings flooded and the floodplains. The key parameters (P) are: the boundary conditions, the number n and the cell size. The following table (Table 3-1 and 3-2) describes the values taken by the key parameters along the sensitivity analysis. 46 Table 3-1: Value of the cell size of the DEMs used to interpolate the LIDAR ground points. 4.8124 4.723 4.622 4.521 4.3920 4.2819 4.1618 4.0417 3.9216 3.815 3.6714 3.5413 3.412 3.2511 3.110 2.949 2.778 2.67 2.46 2.195 1.964 1.853 1.852 1.851 Cell Size (m)Parameter n 6.5845 6.5144 6.4343 6.3642 6.2841 6.240 6.1339 6.0538 5.9737 5.8936 5.835 5.7234 5.6433 5.5532 5.4631 5.3730 5.2829 5.1928 5.127 526 4.9125 5.127 526 4.9125 Cell Size (m)Parameter n 8.1569 8.0968 8.0367 7.9766 7.9165 7.8564 7.7963 7.7362 7.6661 7.660 7.5459 7.4758 7.4157 7.3456 7.2855 7.2154 7.1453 7.0752 7.0151 6.9450 6.8749 6.848 6.7347 6.6546 Cell Size (m)Parameter n 8.2170 Cell Size (m)Parameter n (m 47 Table 3-2: Boundary conditions for the 2-year 5-year and 25-year flood For n equal to one and two, the criteria used to define the cell size of the DEMs could not be respected. Indeed, defining the cell size so that there is in average one ground point per cell would lead to a very small cell size (less than 1.86 meter: 0.98 meter for n=1, 1.39 meter for n=2, and 1.70 for n=3), and the number of points defining a cross section would be too high to be handled with HEC-RAS (see next chapter). Therefore, for the two highest densities of LIDAR ground points (n equal to 1 and 2) the minimum cell size allowed by HEC-RAS was used (1.86 meter). The next chapter provides the theoretical background and the main steps required to construct the hydraulic model. 7.8422.525-year Flood 7.3114.25-year Flood 6.87.82-year Flood Water surface elevation at the downstream cross-section Q (m 3 ) 48 CHAPTER FOUR: CONSTRUCTION OF THE HYDRAULIC FLOOD MODEL 4.1 DESCRIPTION OF THE AREA OF LIDAR SURVEY The LIDAR data were collected in Brownsville, Texas, and in the northern portion of Matamoros, Mexico. These two cities located along the Rio Grande have followed a very different history of urban development. These differences are stamped on the hydrologic landscape. 4.1.1 The hydrologic landscape in Brownsville Brownsville is located 38 km from the mouth of the Rio Grande (Figure 4-1). In addition to the Rio Grande, other major fluvial features in the lower Rio Grande Valley are the Arroyo Colorado and the North Floodway, which fall outside of this study area. The North Floodway was constructed to minimize flood risks from the Rio Grande, and the Arroyo Colorado is a natural drainage feature that is maintained for irrigation purposes. In the Brownsville vicinity, the Brownsville Ship Channel and numerous resacas (remnants of old river channels with considerable linear extent) are the largest hydrologic features. Old oxbow meanders (locally called bancos), levees, irrigation canals, and manmade reservoirs also play a vital hydraulic role. The main hydraulic features are from south to north: the Town Resaca, the North Main Drain, Resaca de la Guerra, and Cameron County Ditch (Figure 4-2). 49 Figure 4-1: Satellite image (SPOT) of Brownsville and Matamoros 50 Figure 4-2: Main resacas and drains in Brownsville with their watersheds Several industrial sites and marsh wetlands are concentrated near the inland terminus of the Brownsville Ship Channel that extends within 3 km of the urban area, northeast of the city. Near to the city, resacas provide multiple benefits, such as collection and storage of local storm runoff, conveyance channels for Rio Grande waters, irrigation and drinking water sources, habitat for wildlife, and recreational attraction. New housing developments are commercializing on aesthetic appeal of resacas-frontage lots. The main resacas in Brownsville include Resaca de la Guerra, Fort Brown Resaca, Resaca del Rancho Viejo, and Town Resaca. The major Bancos in the area (which are smaller than resacas) include Los Tomato Banco, Matamoros Banco and Morales Banco. Rio Grande Town Resaca North Main Resaca de la Guerra Cameron County Ditch 51 The major levees are located along the banks of the Rio Grande and are maintained by the International Boundary and Water Commission. Since the Rio Grande forms the International Boundary with Mexico, the purpose of these levees is to control the flow of the Rio Grande to minimize changes to the stream channel and the boundary demarcation, as well as flood control for the cities. The extensive irrigation canals support year-round agriculture. The source of water includes resacas, man-made storage reservoirs, and the Rio Grande. Sixteen major pumping stations lift water from the Rio Grande into conveyance canals for the irrigation districts. With minimum topographic relief in the area, many of these hydrologic features flow in both directions depending on the action of pumping stations. 4.1.2 The hydrologic landscape in Matamoros The natural landscape of Matamoros, south of the Rio Grande, was once characterized by an abundance of abandoned river meanders, called resacas. As in Brownsville, the resacas near the river are the most recent. Historically the resacas played an important role in settling the area as by suppling water for agriculture and providing transportation routes. Modern Matamoros has covered up with urbanization all but three of these resacas. The city?s water department, Junta de Aquas y Drenajes, realizes the importance of resacas, along with smaller bancos, for flood control and aesthetic appeal; however, their loss to urbanization cannot easily be reversed. Over the two last centuries the growth of population forced the paving over of the resacas: San Pablo, del Bravo, Lui Cabaleros, and Laguna de San Fransisco. Several intermittent lagoons, oxbows, and streams still remain on the outskirts of Matamoros. Resaca de las Rusias is northwest of Matamoros, Arroyo del Tigre is south of the city, and the lagoons of Llano Salado and Palangana lie to the southeast. During most of the year these features are dry and only fill with water during rainy periods. Since 52 these features are on the coastal plain and subject to high evaporation rates, their water is saline. Other important drainage features in Matamoros area include agricultural and open sewage drainage canals. The main drainage canals are the Dren Viente de Noviembre (east), Dren E-32 Izquierdo (west), Dren Emisor Pluvial (west), Dren Principal, Dren de las Vacas, and Dren del Mar (southwest), all of which transport water from the north to the south. Unfortunately, due to the lack of environmental regulations to maintain these canals, the discharge capacity of these canals during floods may limit their effectiveness. The most obvious difference in the surfacial hydrography of the two cities is that Brownsville has left many resacas in their natural form, whereas in Matamoros dense urbanization has filled in the abandoned channels. Also an extensive floodway system has been constructed for both irrigation and runoff control on the US side, whereas in Mexico many canals are used for sewage and trash disposal. 4.1.3 Climate and historical floods The study area lies in a sub-tropical humid climate that is strongly influenced by Gulf-related weather activity. The average annual rainfall for the area is 665 mm, with two seasons of maximum precipitation. Precipitation from May to June represents 21% of the annual total and from August to October another 42%. During the winter evapotranspiration may exceed precipitation. Late summer hurricanes and tropical storms initiate much of the area?s precipitation. Therefore, it is common to see heavy rainfall amounts of 250 to 500 mm in a short time period (Mendoza 2001). Since the area lies in a low relief landscape, only 5 to 6 meters above sea level, with a groundwater table only one meter from the surface, short duration rainfall amounts of 50 mm or more cause flood damage in both the Brownsville and Matamoros urban areas. Here are some examples: 53 ? In 1967, Hurricane Beulah flooded 33% of the Brownsville-Harlingen area and dropped over 300 mm of rain in Matamoros. The majority of Matamoros was underwater, including the main electrical substation, cutting off outside communication. ? In September 1980, Hurricane Allen struck north of Brownsville and caused coastal and inland flooding from storm surges, heavy rainfall, and wind damage. ? In September 1984, in four days, 467 mm of rain fell, causing widespread flooding. ? In 2001, over 70% of Matamoros was flooded heavily. 4.1.4 The North Main Drain The study area is a portion of the North Main Drain (Figure 4-3) located in Brownsville (Texas). A small portion of the drain located between two junctions was chosen for this study. At least three reasons can be given to justify this choice: ? The North Main Drain is a fairly conventional channel system. The area of resacas was not chosen because the hydraulics of a resaca is much more complicated. Indeed, resacas are a succession of pools connected by pipes and pumps. ? A previous flood study on the North Main Drain was conducted by Rice University, with a 30 meter cell size DEM. It could be interesting to compare the floodplains obtained by Rice University and the one obtained by this study (Bedient, 2003). ? The study area is located between two junctions because flow data were not available for the reaches beyond these junctions. 54 ? Only a small portion of the North Main Drain was chosen to reduce the computation time since a large number of simulations were run for the sensitivity analysis. Figure 4-3: North Main Drain near Boca Chica Boulevard 4.2 ONE DIMENSIONAL HYDRAULIC MODELING THEORY This section gives a short overview of the theoretical concepts used by a one- dimensional hydraulic model. The principles used to compute the water surface profile, as well as the basic equations used are given. 55 4.2.1 Energy equation The model developed for this study is a one dimensional flow model, which means that a sequence of cross sections of the river channel were selected and the flow of the river at these cross sections is assumed to be perpendicular to these cross sections. The cross sections of the North Main Drain were extracted from a DEM and the water surface profiles are computed from one cross section to the next one by solving the Energy equation with an iterative procedure. The Energy equation reads: Where: 21 , YY = depth of water at cross section 21 , ZZ = elevation of the main channel inverts 21 , VV = average velocities (total discharge/total flow area) 21 ,?? = velocity weighting coefficients g = gravitational acceleration e h = energy head loss The energy head loss ( e h ) between two cross sections is comprised of friction losses and contraction and expansion losses. The equation of the energy head loss is as follows: Where: L = discharge weighted reach length f S = representative friction slope between two sections C = expansion and contraction loss coefficient The distance weighted reach length is calculated as: e h g V ZY g V ZY +++=++ 2 2 2 11 11 2 22 22 ?? (3) g V g V CSLh fe 22 2 11 2 22 ?? ?+= (4) 56 Where: robchlob LLL ,, = cross section reach length specified for flow in the length overbank, main channel and right overbank, respectively robchlob QQQ ,, = arithmetic average of the flows between sections for the left overbank, channel and right overbank, respectively The determination of the total conveyance and the velocity coefficient for cross section requires that flow be subdivided into units for which the velocity is uniformly distributed. The approach used in HEC-RAS is to subdivide flow in the over bank area using the input cross section n-value break points (location where n value changes) as the basis for division (Figure 4-4). Conveyance is calculated within each subdivision from the following form of Manning?s equation (based on English units): Where K = Conveyance for subdivision n = Manning?s roughness coefficient for subdivision A = flow area subdivision R = hydraulic radius for subdivision (area/ wetted perimeter) Then the HEC-RAS program sums the incremental conveyances in the overbank areas to obtain a conveyance for the left overbank and the right overbank (HEC RAS Hydraulic References). The main channel conveyance is normally computed using a single conveyance element. The total conveyance for the cross section is obtained by summing the three conveyances (Left, channel and right). robchlob robrobchchloblob QQQ QLQLQL L ++ ++ = (5) 2/1 f SKQ = (6) 3/2 486.1 RA n K = (7) 57 Figure 4-4: Cross section of the North Main Drain subdivided into units with different n values 4.2.2 Boundary conditions Boundary conditions are necessary to establish the initial water surface level at the end of the river system (upstream and downstream). For this study, a subcritical Station (meters) Elevation (meters) n value 58 steady state flow regime was used. In that case, boundary conditions are necessary only at the downstream end of the river. For a steady state regime, the boundary condition consists of the flow and height of water at the downstream cross section of the channel. The above equations are solved for each cross section, successively from upstream to downstream, in order to determine the flow and the height of water at each cross section. 4.2 HEC-RAS HYDRAULIC MODEL This section describes how the model was set up in HEC-RAS. The objective of the simulations is to compute water surface elevations at all locations of interest for a given flow (steady flow simulation). The data needed to perform these simulations are divided into the following categories: ? Geometric data ? Steady flow data 4.2.1 HEC-RAS Geometric data Geometric data include a river system schematic and cross section data. The geometric data to be defined to complete the creation of the geometric file are: ? Stream centerlines, which have an associated, attribute table that defines river and reach identifiers. ? Stream banks: left and right banks of the streams defined by lines parallel to the stream centerlines. ? Left bank is defined as the left line when looking downstream. ? Overbanks flowpath: used to derive downstream reach length. ? Land use to derive the roughness value of the terrain. ? Cross section lines 59 No Hydraulic or geometric data concerning bridges and culverts were available for this study. But the idea is to get representative floodplains, not precisely accurate ones. The river system schematic was extracted from the DEM with the highest resolution (1 meter cell size) using the HEC-GeoRAS 3.1.1 preprocessor and then reused for the other simulations with lower resolution DEMs. Indeed, within the range of resolution considered for the sensitivity analysis (1.85 meter cell size to 8.21 meter cell size), the definition of Stream centerlines, Stream banks and Overbanks flowpath were not modified by resolution changes. However, the changes in resolution of the DEMs are critical for cross section extraction. Therefore, cross sections were extracted automatically from each DEM using the Model Builder. 4.2.1.1 The river system schematic Although HEC-RAS can route through a dendritic network of reaches, for this study, a single reach was defined with its banks and flowpath (Figures 4-4 and 4-5). The following figures summarize the characteristics of these feature classes. 2 Kilometers 2.2 Kilometers 60 Figure 4-5: The river system schematic (stream banks, stream flowpath and stream centerline) and their table of contents. 4.2.1.2 The cross sections Boundary geometry for the analysis of one-dimensional flow in natural streams is specified in terms of ground surface profiles (cross sections) and the measured distance between them (reach lengths). Cross sections are located at intervals along a stream to characterize the flow carrying capability of the stream and its adjacent floodplain. Cross- sections extend across the entire floodplain and are perpendicular to the anticipated flow lines. A total of 48 cross sections were defined for the model (Figure 4-6). 61 Figure 4-6: A total of 48 cross sections were used to compute the flow of the North Main Drain. The green area on the east of the North Main Drain is a Resaca which is not included in the model. The 2-year, 5-year and 25-year floodplain do not overlap the Resaca. Each data point on the cross sections is given a station number corresponding to the horizontal distance from a starting point on the left of the cross section when looking downstream. Up to 500 data points can be used to describe each cross section. Problems due to this limitation appear when a high resolution DEM is used to extract the 3D cross sections. Indeed, the methodology described in the previous chapter stipulates that the cell size of the DEM be chosen so that there is on average one LIDAR ground point per cell. The highest densities of LIDAR ground point (104 points/100m 2 and 52 Resaca 62 points/100m 2 ) would yield to a 1.01 meter cell size grid and a 1.56 cell size grid, which would give more than 500 data points on the longest cross section of the model. Therefore, the LIDAR ground points datasets with a density equal to 104 points/100m 2 and 52 points/100m 2 were interpolated with the smallest cell size allowed by HEC-RAS (1.86 meter) which gives 500 data points on the longest cross section. It would have been interesting to study the 100-year flood event because the Federal Insurance Management Agency uses the 100-year flood to create the National Flood Insurance Rate Maps. Due to the limitation of the HEC-RAS software (a maximum of 500 points to describe a cross section), the 100-year flood event was not studied although data were available. Indeed, the cross sections length could not cover the 100- year floodplain extent and be described by less than 500 points for the highest densities of LIDAR ground points. There is an option in HEC-RAS to filter ?smartly? the points on a cross-section but no good results could be obtained. The friction loss is estimated using Manning?s n value. The selection of an appropriate n value is very significant to the appropriate water surface profiles. The value of Manning?s is highly variable and depends on a number of factors including: surface roughness including vegetation characteristics, channel irregularities, structural obstructions, size and shape of the channel. Manning?s n coefficients were calibrated using the observed surface profile. Fifteen types of landuse classes were used (Figures 4- 7 and 4.8). 63 Figure 4-7: Shapefile of landuse with the cross sections of the North Main Drain 64 Figure 4-8: Table of contents of the laduse shapefile used to extract the Manning's n values An extensive compilation of the Manning?s n values for typical channel and floodplains was found in ?Open Channel Hydraulics? (Chow 1959) and is displayed in Table 4.1. 65 Table 4.1: Manning's n values (Chow 1959) 4.2.2 Calibration of the model The model was calibrated for the 2-year, 5-year and 25-year flood event. As has already been mentioned, one of the options available in HEC-RAS to define the boundary conditions for a steady state simulation is to fix the flow and the height of water at the downstream cross section of the channel. This option was used to calibrate the model. 66 The data available for the calibration process were taken from the Flood protection Plan (November 1996, revised in 1997) and describe the water surface profiles along the North Main Drain of the 2-year, 5-year and 25-year flood event. The objective of the calibration is to find the flow and the height of water at the downstream cross section for the 2-year, 5-year and 25-year flood so that the computed water surface profiles fit the water surface profiles produced in this study. The two main steps of the calibration are: ? First, for a given flood event, the height at the downstream cross section of the channel was read on the water surface profile of the North Main Drain given by the study. ? Second, a flow at the downstream cross-section of the channel was determined by trial and error so that the computed water surface profile of the North Main Drain fit the displayed water surface profile (Figure 4-9). 67 Figure 4-9: Calibration of the hydraulic flood model by adjusting the flow at the downstream cross section so that the computed water surface profile fits with the displayed water surface profile. The boundary conditions obtained with the calibration process are summarized in Table 4-2. Table 3-2: Boundary conditions (flow and water surface elevation at the downstream cross section of the North Main Drain) for the 2-year, 5-year and 25-year flood 0 500 1000 1500 2000 2500 6.3 6.4 6.5 6.6 6.7 6.8 6.9 7 7.1 7.2 7.3 7.4 7.5 7.6 7.7 7.8 7.9 8 Distance along the channel from downstream in meters Wa t e r s u r f a c e el ev a t i o n i n m e t e rs 2-year flood water surface elevation simulation 5-year flood water surface elevation simulation 25-year flood water surface elevation simulation 2-year flood water surface elevation 5-year flood water surface elevation 25-year flood water surface elevation 7.8422.5025-year flood 7.3114.205-year flood 6.907.802-year flood H (m)Q (m3/s) 68 Once calibrated, the HEC-RAS model was used to study its sensitivity to the density of LIDAR ground points. The results are analyzed and discussed in the next chapter. 69 CHAPTER FIVE: RESULTS AND DISCUSSION Three types of results were analyzed: the influence of the density of LIDAR ground points on the extent and the location of the computed floodplains and watersheds and the influence of the density of the LIDAR ground points on the number of buildings flooded. It was also interesting to look at the evolution of the computation time of the simulation as the density of LIDAR data decreases. 5.1 COMPARISON OF FLOODPLAIN AND WATERSHED DELINEATION AS THE DENSITY OF LIDAR POINTS DECREASES Comparison of floodplains or watersheds involves the same principle since both are represented as polygons. The sets of polygons (floodplains or watersheds) were compared with the polygon obtained with the highest density of LIDAR points (which is considered to be the most accurate). To compare these polygons and to take into account the differences in extent and location, the error was defined as described in Chapter Two. An example is repeated below (Equation 1 and Figure 5-1): Figure 5-1: Illustration of the definition of the error 100 )( )(2)()( (%) ? ???+ = PolygonREFArea PolygonREFPolygonAreaPolygonREFAreaPolygonArea Error (2) Polygon REF Polygon ErrorIntersection 70 5.1.1 Floodplain results For a given flood event, the flood polygons obtained with the 70 different densities of LIDAR ground points were compared with the flood polygon obtained using the highest density of LIDAR ground points. An example of the 5-year floodplain computed with the highest density of LIDAR ground points and the lowest density of LIDAR ground points is given on Figure 5-2. More results are presented in appendices. Figure 5-2: Example of the differences in location and extent between the 25-year floodplain computed with highest (left) and lowest (right) density of ground LIDAR points The following graphs (Figures 5-3, 5-4 and 5-5) represent for the 2-year, 5 year- flood and 25-year flood events, and the evolution of the error as the density of LIDAR points changes. 71 Figure 5-3: Calculation of the error in the delineation of the flood polygons as the density of LIDAR points changes for the 2-year flood 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 0 50 100 150 200 250 300 350 400 450 500 Density of LIDAR points in points/m2 E rro r i n % 72 Figure 5-4: Calculation of the error in the delineation of the flood polygons as the density of LIDAR points changes for the 5-year flood 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 Density of LIDAR points in points/m2 E r ro r i n % 73 Figure 5-5: Calculation of the error in the delineation of the flood polygons as the density of LIDAR points changes for the 25-year flood For the 25-year flood, only the first 26 simulations were completed. Indeed, as the resolution of LIDAR data decreases, the extent of the computed floodplain increases. In the case of the North Main Drain, the computed extent became equal to the maximum extent defined by the cross sections and was then limited by the width of the cross sections. For the three flood events, there is a break point around 0.05 points/m 2 that corresponds to a little more than 5 points per 100m 2 . This break point is due to the small width of the channel. Indeed, the average distance between points corresponding to this density is 4.5 m, which is very close to the width of the North Main Drain. Therefore, for 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0.18 0.2 Density of LIDAR points in points/m2 E r ro r i n % 74 this density and lower densities of points, flow obstructions appear on the 3D cross- sections derived from LIDAR (Figure 5-6), and the inundation polygons expand to unrealistic proportions (characterized by a high error). Figure 5-6: Apparition of flow obstructions for low densities of LIDAR points The appearance of flow obstructions on the 3D cross sections are characterized by a high variation of water surface elevation or a discontinuity of the computed water surface profile (Figure 5-7). 0 50 100 150 200 250 300 350 6.5 7 7.5 8 8.5 9 9.5 10 10.5 11 Distance in meter from the left to the right of the cross section when looking downstream El e v a t i o n i n m e t e r 3D cross section obtained with a density of LIDAR ground points of 1 point/m2 3D cross section obtained with a density of LIDAR ground points of 1point/100m2 25 meters 75 Figure 5-7: The appearance of flow obstructions on the 3D cross sections when using a low density of LIDAR ground points creates discontinuities in the computed water surface profiles due to artificial flow obstruction For a density lower than 0.055 point/m 2 the LIDAR data cannot be used to delineate the 2-year or 5-year flood event. At densities greater than 5 ground points per 100 m 2 , the influence of LIDAR point density is weak, and the accuracy of the floodplain 0 500 1000 1500 2000 2500 6 7 8 9 10 11 Diss8 Plan: Plan 01 Main Channel Distance (m) E l ev at ion (m) Legend EG PF 1 WS PF 1 Crit PF 1 Ground 1 0 500 1000 1500 2000 2500 6 7 8 9 10 11 Diss10 Plan: Plan 01 Main Channel Distance (ft) Elevation (ft) Le EG WS Cri G 1 0 500 1000 1500 2000 2500 6.0 6.5 7.0 7.5 8.0 8.5 9.0 9.5 10.0 Diss5 Plan: Plan 01 Main Channel Distance (m) Elevation (m) Legend EG PF 1 WS PF 1 Crit PF 1 Ground 1 0 500 1000 1500 2000 2500 6 7 8 9 10 11 Diss6 Plan: Plan 01 Main Channel Distance (m) E l evat i on ( m ) L E W C G 1 0 500 1000 1500 2000 2500 6.0 6.5 7.0 7.5 8.0 8.5 9.0 9.5 10.0 Diss3 Plan: Plan 01 Main Channel Distance (m) Elevat ion (m ) Legend EG PF 1 WS PF 1 Crit PF 1 Ground 1 0 500 1000 1500 2000 2500 6.0 6.5 7.0 7.5 8.0 8.5 9.0 9.5 10.0 Diss4 Plan: Plan 01 Main Channel Distance (m) E l ev at ion (m) Lege EG P WS P Crit P Grou 1 El ev at i on (m) Elevation (ft) El ev at i on (m) El ev at i on (m) El ev at i on (m) Elevation (ft)Elevation (ft)Elevation (ft) Elevation (m) E l evat i on ( m ) Elevat ion (m ) E l ev at ion (m) Elevation (m) E l evat i on ( m ) Elevation (m)Elevation (m)Elevation (m) E l evat i on ( m ) E l evat i on ( m ) Elevat ion (m ) E l ev at ion (m) Elevat ion (m ) E l ev at ion (m) 1.04 point/m 2 0.10 point/m 2 0.065 point/m 2 0.055 point/m 2 0.025 point/m 2 0.015 point/m 2 76 resulting from a simulation depends upon the initial conditions. For example, for a density of 0.1 point/m 2 , there is almost a 50 percent error for the 2-year flood event, whereas there is a 0.2 percent error for the 5-year flood event and only a 0.1 percent error for the 25-year flood. 5.1.2 Results concerning watersheds No clear trend was observed for the delineation of watersheds as the density of LIDAR points decreases (Figure 5-8). This failure to determine a relationship between ground point density and accurate watershed delineation could be explained by the fact that the area of study is extremely flat. Figure 5-8: Calculation of the error in the delineation of watersheds as the density of LIDAR points changes 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 Density of LIDAR points in points/m2 E r ro r i n % 77 Figure 5-9 shows an example of the watersheds computed with the highest density of LIDAR ground point and the lowest density of LIDAR ground points. More results are presented in appendices. Figure 5-9: Watershed of the North Main Drain computed with the highest density of LIDAR ground points (104 points per 100m 2 , left) and the lowest density of LIDAR ground points (1.5 points per 100m 2 , right) 5.2 COMPARISON OF THE NUMBER OF BUILDINGS FLOODED AS THE DENSITY OF LIDAR POINTS DECREASES The number of buildings flooded using the highest density of LIDAR data was taken as a reference. Two types of situations can occur as the density of LIDAR ground points decreases: ? A building can be ?mistakenly? flooded, which means that the building is not flooded when using the highest density of LIDAR ground points 78 (which theoretically describes the most accurate case) but it is flooded when using a lower density. ? A building can be ?mistakenly? not flooded, which means that the building is flooded when using the highest density of LIDAR ground points but it is not flooded when using a lower density. Both cases were analyzed. It turns out that the number of buildings ?mistakenly not flooded is always zero as the density of LIDAR points increases. The following graphs (Figures 5-10, 5-11, 5-12) show calculations of the number of buildings mistakenly flooded as the density of LIDAR ground points decreases. Figure 5-10: Calculation of the number of buildings mistakenly flooded by the 2-year flood as the density of LIDAR ground points changes 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 0 100 200 300 400 500 600 700 800 Density of LIDAR points in points/m2 N u m b e r of bu i l di ng s m i s t ak e n l y f l o oded 79 Figure 5-11: Calculation of the number of buildings mistakenly flooded by the 5-year flood as the density of LIDAR ground points changes 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 0 50 100 150 200 250 300 350 400 450 500 550 600 650 700 Density of LIDAR points in points/m2 N u m b e r of bu i l di ng s m i s t ak e n l y fl o oded 80 Figure 5-12: Calculation of the number of buildings mistakenly flooded by the 25-year flood as the density of LIDAR ground points changes The shape of the curve reflects the same trend as the previous graphs, but it is significant to note that a too low density of LIDAR ground points can lead to a very high overestimation of the number of buildings flooded. 5.3 COMPUTING RESOURCES It is interesting to look at the evolution of the processing time as the density of LIDAR points decreases because computation time is a key cost factor besides the cost of data collection and preparation. The sensitivity analysis was run on an Intel? Pentium? 4, CPU 2.53 GHz with 1 MB of RAM. 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 0 10 20 30 40 50 60 70 80 90 100 Density of LIDAR points in points/m2 N u m ber o f bui l d i n g s m i s t ak enl y f l ooded 81 The computation time of some of the loops were recorded and plotted on the following graph (Figure 5-13): Figure 5-13: Evolution of the computation time required to decimate the LIDAR ground points and compute the floodplain as density of LIDAR ground points decreases. The computation time varies considerably (427%) with the density of LIDAR ground points. Therefore it is an important parameter to take into account. It took 3 days and 11 hours to run 70 flood simulations for a given flood event. 70 simulations were run for the 2-year and 5-year flood event and 26 simulations were run for the 25-year flood event. Besides flood simulations, the watersheds were delineated for the seventy different densities of LIDAR ground points. It took a bit more than 10 days to run the entire sensitivity analysis, which largely justifies the need to automate the process. 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 0 20 40 60 80 100 120 140 160 180 200 220 Density of LIDAR ground points on points/m2 C o m p u t at i on t i m e i n m i nut es 82 CHAPTER SIX: CONCLUSIONS AND RECOMMENDATIONS 6.1 CONCLUSIONS The density of LIDAR ground points is a crucial parameter for the purposes of flood modeling. Results show that for the study area, the density of LIDAR ground points could be reduced to 0.055 point/m 2 (which corresponds to a bit more than 5 points per 100 m 2 or a post-spacing of 4 meters) without impacting the accuracy of the results of the flood model (only a few percent error). Below the threshold of 0.055 point/m 2 the sensitivity of the model to the density of LIDAR ground points depends on the initial condition used for the flood simulation. At lower point densities, flow obstructions appear on the 3D cross-sections derived from LIDAR, and the inundation polygons expand to unrealistic proportions. It is interesting to see that we could reduce the density of the LIDAR data (and the cost) to 0.055 point/m 2 and obtain the same accuracy with respect to the flood hazard simulations. This density of 0.055 point/m 2 is linked to the width of the river. Indeed, a density of LIDAR ground points of 0.055 point/m 2 corresponds to 0.23 point/m (or a 4.34 meter post-spacing), and the North Main Drain?s width is approximately 25 meters (from the left bank to the right bank), which means that a at least 4.75 points (or let?s say 5 points) in the channel are required to appropriate describe a cross section. This result allows transferring the results of this study (which was done on a drain) to a large river. Concerning watershed delineation, no trend was observed as the density of LIDAR points decreases. 83 6.2 RECOMMENDATIONS Brownsville Texas is located in a very flat area near the mouth of the Rio Grande. Conclusions drawn from the study may be applicable to flood hazard mapping in many other regions having low surface relief. A similar study could also be done in a region with higher surface relief and then compared with the current study. Because it appears that the width of the channel is a key parameter to estimate the value of the density of LIDAR ground points at the break point, another study with a larger river such as the Rio Grande could be done to verify this hypothesis. The study was performed with a one-dimensional flood model. Flood impact assessment has traditionally relied on the use of one-dimensional models such as HEC- RAS. These models are able to accurately predict flood levels and discharges in applications where the basic assumptions of one-dimensional flow remain valid. Difficulties can arise, however, in urban areas such as the city of Brownsville, where the flow paths on the floodplain can become ill-defined due to the incidence of man made structures such as houses, roadways, bridges, embankments and levee banks, etc. In particular, problems arise where the analysis is required at a level of detail involving flows around individual buildings and structures where a one dimensional model is limited in its capacity to accurately represent the resulting complex flow paths. The aim of 2D modeling is not normally to model extensive reaches of a river system but, instead, to focus on specific segments of river and floodplain (such as the portion of the North Main Drain studied in this project), where detailed analysis of the flow conditions is required. Therefore, in order to make full use of the possibilities offered by high- resolution topographical data, such as that obtained by LIDAR, a two dimensional flood model could be developed, and a similar study conducted. Some simulations could be made with the buildings since the two dimensional flow can be computed around the 84 buildings. The advantage of two dimensional flood models is now largely recognized (McCowan, 1999). In this project, LIDAR ground points were interpolated into a grid to describe the terrain surface because ArcGIS cannot handle the volume of TIN data involved. The advantage of grid interpolations is that complex surface model algorithms are better suited for a continuous array of evenly spaced data rather than an irregular spatial distribution. It would have also been interesting, for a given density of LIDAR ground points, to compare the floodplain computed with a grid and the floodplain computed with a TIN (Triangulated Irregular Network). Two advantages with TINs are that it does not interpolate beyond the distribution of the data, and the map is forced to fit the control points. TINs are also well suited for the purpose of generating a more highly resolved structure on which to drape or hang feature layers. Additional tests could be done through collection or processing of LIDAR data. For example, the sensitivity test could be rub on the "last-return" data would be more representative of actual data collection made from an aircraft than the analysis done on the already filtered ground point data. These DEMs would be better predictors of the actual performance of an aerial LIDAR sensor at various altitudes than the DEMs created from pre-filtered ground points. It would also be interesting to compare actual LIDAR datasets collected at different altitude. 85 APPENDIX A: REULTS A.1 WATERSHEDS DELINEATION RESULTS (TO SAVE SPACE, ONLY A PORTION OF THE RESULTS ARE PRESENTED) n=1 1.0397 point/m 2 n=10 0.1039 point/m 2 n=20 0.0519 point/m 2 n=30 0.0346 point/m 2 n=40 0.0207 point/m 2 n=50 0.0259 point/m 2 86 n=60 0.0173 point/m 2 n=70 0.0148 point/m 2 87 A.2 TWO-YEAR FLOOD RESULTS (TO SAVE SPACE, ONLY A PORTION OF THE RESULTS ARE PRESENTED) n=1 1.039 point/m 2 n=10 0.1039 point/m 2 n=20 0.051 point/m 2 n=30 0.034 point/m 2 88 n=40 0.025 point/m 2 n=50 0.020 point/m 2 n=60 0.017 point/m 2 n=70 0.014 point/m 2 89 A.3 FIVE-YEAR FLOOD RESULTS (TO SAVE SPACE, ONLY A PORTION OF THE RESULTS ARE PRESENTED) n=1 1.039 point/m 2 n=10 0.1039 point/m 2 n=20 0.051 point/m 2 n=30 0.034 point/m 2 90 n=40 0.025 point/m 2 n=50 0.020 point/m 2 n=60 0.017 point/m 2 n=70 0.014 point/m 2 91 A.4 TWENTY FIVE-YEAR FLOOD RESULTS (TO SAVE SPACE, ONLY A PORTION OF THE RESULTS ARE PRESENTED) n=1 1.039 point/m 2 n=20 0.051 point/m 2 n=10 0.1039 point/m 2 n=26 0.039 point/m 2 92 APPENDIX B: CODE B.1 DECIMATE TOOL CODE Private Function DisseminateLinear(dbpath As String, N As Long) As Boolean Dim dbs As Database Dim strSQL As String Dim i As Integer DisseminateLinear = False Set dbs = DBEngine.Workspaces(0).OpenDatabase(dbpath) For i = N To N strSQL = "SELECT LIDAR_without_bridge.ID, LIDAR_without_bridge.X, LIDAR_without_bridge.Y, LIDAR_without_bridge.Z " strSQL = strSQL & " INTO tbl" & i & " from LIDAR_without_bridge " strSQL = strSQL & " WHERE ([ID] Mod " & i & ")=0" Debug.Print (strSQL) dbs.Execute (strSQL) Next i DisseminateLinear = True End Function Private Function DisseminateSquare(dbpath As String, N As Long) As Boolean Dim dbs As Database Dim strSQL As String Dim i As Integer DisseminateSquare = False Set dbs = DBEngine.Workspaces(0).OpenDatabase(dbpath) For i = 1 To N strSQL = "SELECT LIDAR_without_bridge.ID, LIDAR_without_bridge.X, LIDAR_without_bridge.Y, LIDAR_without_bridge.Z " strSQL = strSQL & " INTO tblsq" & i ^ 2 & " from LIDAR_without_bridge " strSQL = strSQL & " WHERE ([ID] Mod " & i ^ 2 & ")=0" Debug.Print (strSQL) 93 dbs.Execute (strSQL) Next i DisseminateSquare = True End Function Public Function AccessDLL(Thedbpath, TheNvalue, TheCriteria) As Boolean 'To make connectivity with script tool, pass script-parameters and call main function AccessDLL = False 'If successful it becomes True 'Set global variables equal to the input variables Dim sThedbpath As String Dim lNvalue As Long Dim sTheCriteria As String sThedbpath = Thedbpath 'From script-variant to string lNvalue = TheNvalue 'From script-variant to string sTheCriteria = TheCriteria 'From script-variant to string 'Call the main function to continue processing If sTheCriteria = "Linear" Then If DisseminateLinear(sThedbpath, lNvalue) Then 'If main function succeeds it becomes True AccessDLL = True 'Indicates success End If ElseIf sTheCriteria = "Quadratic" Then If DisseminateSquare(sThedbpath, lNvalue) Then 'If main function succeeds it becomes True AccessDLL = True 'Indicates success End If End If End Function B.2 RUN HEC-RAS TOOL CODE Public sTheRASprj As String Public The.SDF As String Public sTheGeo As String Private HecRas As RAS.HECRASController 94 Private Function CallRAS() As Boolean CallRAS = False 'If code is NOT sucessful the script will know Dim strArgs As String Dim sStartTime As Single strArgs = sTheRASprj strArgs = Trim(strArgs) Set HecRas = New RAS.HECRASController HecRas.ProjectName = strArgs HecRas.RunSteady sStartTime = Timer Do While ((Timer - sStartTime) < 4) DoEvents If ((Timer - sStartTime) > 10) Then Exit Do Loop HecRas.ExportGIS '''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''' 'Check for existence of new .sdf file '''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''' 'Create new file system object 'Dim fso As FileSystemObject 'Set fso = New FileSystemObject 'Build the string for the new filename Dim .SDFFile As String 'The .SDFfile is the same Project filename with a .DSC extension .SDFFile = The.SDF 'Stop On Error Resume Next HecRas.Quit CallRAS = True 'Thus, If code is sucessful the script will know End Function 95 Public Function AccessPierre(TheRASprj, The.SDF, TheGeo) As Boolean 'To make connectivity with script tool, pass script-parameters and call main function AccessPierre = False 'Set global variables equal to the input variables sTheRASprj = TheRASprj 'From script-variant to string sTheGeo = TheGeo 'Call the main function to continue processing If CallRAS Then AccessPierre = True 'Indicates success End If End Function 'Create the Geoprocessor object Set gp = WScript.CreateObject("esrigeoprocessing.GPDispatch.1") 'Load required toolboxes... gp.AddToolbox "C:/Model_Builder_Lidar/SensitivityAnalysisTool.tbx" 'Local variable... lidarpointsamy_mdb = "C:\Model_Builder_Lidar\Input_Sensitivity_analysis\lidarpointsamy.mdb" dim Output_Watershed(70) for i=1 to 70 Output_Watershed(i) = "C:/Model_Builder_Lidar/Output_Sensitivity_Simulation/Watersheds/watershed"&i next dim Output_watershedpoly(70) for i=1 to 70 Output_watershedpoly(i) = "C:/Model_Builder_Lidar/Output_Sensitivity_Simulation/Watersheds/watershedpoly"&i next dim Intersect_watershedpoly(70) 96 for i=1 to 70 Intersect_watershedpoly(i) = "C:/Model_Builder_Lidar/Output_Sensitivity_Simulation/Watersheds/interwatershedpol y"&i next dim waterchrono(70) for i=1 to 70 waterchrono(i) = "Z:\outputmodelbuilder_chrono\Watersheds\watershedpoly"&i&".shp" next B.3 LOOPING TOOL CODE dim output_XS(70) for i=1 to 70 output_XS(i) = "C:/Model_Builder_Lidar/Output_Sensitivity_Simulation/3D_feature_for_HecRas/XSect ion3D"&i next dim Shapefile(70) for i=1 to 70 Shapefile(i) = "C:/Model_Builder_Lidar/Output_Sensitivity_Simulation/LidarPoints/LidarPoints"&i next dim LidarRaster(70) for i=1 to 70 LidarRaster(i) = "C:/Model_Builder_Lidar/Output_Sensitivity_Simulation/LidarRaster/raster"&i next dim iteration(70) for i=1 to 70 iteration(i)=i next dim outputcellsize2(70) outputcellsize2(1)=1.85 outputcellsize2(2)=1.85 outputcellsize2(3)=1.85 outputcellsize2(4)=1.96 97 outputcellsize2(5)=2.19 outputcellsize2(6)=2.40 outputcellsize2(7)=2.60 outputcellsize2(8)=2.77 outputcellsize2(9)=2.94 outputcellsize2(10)=3.10 outputcellsize2(11)=3.25 outputcellsize2(12)=3.40 outputcellsize2(13)=3.54 outputcellsize2(14)=3.67 outputcellsize2(15)=3.80 outputcellsize2(16)=3.92 outputcellsize2(17)=4.04 outputcellsize2(18)=4.16 outputcellsize2(19)=4.28 outputcellsize2(20)=4.39 outputcellsize2(21)=4.50 outputcellsize2(22)=4.60 outputcellsize2(23)=4.70 outputcellsize2(24)=4.81 outputcellsize2(25)=4.91 outputcellsize2(26)=5.00 outputcellsize2(27)=5.10 outputcellsize2(28)=5.19 outputcellsize2(29)=5.28 outputcellsize2(30)=5.37 outputcellsize2(31)=5.46 outputcellsize2(32)=5.55 outputcellsize2(33)=5.64 outputcellsize2(34)=5.72 outputcellsize2(35)=5.80 outputcellsize2(36)=5.89 outputcellsize2(37)=5.97 outputcellsize2(38)=6.05 outputcellsize2(39)=6.13 outputcellsize2(40)=6.20 outputcellsize2(41)=6.28 outputcellsize2(42)=6.36 outputcellsize2(43)=6.43 outputcellsize2(44)=6.51 outputcellsize2(45)=6.58 outputcellsize2(46)=6.65 outputcellsize2(47)=6.73 outputcellsize2(48)=6.80 outputcellsize2(49)=6.87 98 outputcellsize2(50)=6.94 outputcellsize2(51)=7.01 outputcellsize2(52)=7.07 outputcellsize2(53)=7.14 outputcellsize2(54)=7.21 outputcellsize2(55)=7.28 outputcellsize2(56)=7.34 outputcellsize2(57)=7.41 outputcellsize2(58)=7.47 outputcellsize2(59)=7.54 outputcellsize2(60)=7.60 outputcellsize2(61)=7.66 outputcellsize2(62)=7.73 outputcellsize2(63)=7.79 outputcellsize2(64)=7.85 outputcellsize2(65)=7.91 outputcellsize2(66)=7.97 outputcellsize2(67)=8.03 outputcellsize2(68)=8.09 outputcellsize2(69)=8.15 outputcellsize2(70)=8.21 dim outputcellsize(70) outputcellsize(1)=6.07 outputcellsize(2)=6.07 outputcellsize(3)=6.07 outputcellsize(4)=6.44 outputcellsize(5)=7.20 outputcellsize(6)=7.88 outputcellsize(7)=8.51 outputcellsize(8)=9.10 outputcellsize(9)=9.65 outputcellsize(10)=10.18 outputcellsize(11)=10.67 outputcellsize(12)=11.15 outputcellsize(13)=11.60 outputcellsize(14)=12.04 outputcellsize(15)=12.46 outputcellsize(16)=12.87 outputcellsize(17)=13.27 outputcellsize(18)=13.65 outputcellsize(19)=14.03 outputcellsize(20)=14.39 outputcellsize(21)=14.75 outputcellsize(22)=15.10 99 outputcellsize(23)=15.44 outputcellsize(24)=15.77 outputcellsize(25)=16.09 outputcellsize(26)=16.41 outputcellsize(27)=16.72 outputcellsize(28)=17.03 outputcellsize(29)=17.33 outputcellsize(30)=17.63 outputcellsize(31)=17.92 outputcellsize(32)=18.21 outputcellsize(33)=18.49 outputcellsize(34)=18.77 outputcellsize(35)=19.04 outputcellsize(36)=19.31 outputcellsize(37)=19.58 outputcellsize(38)=19.84 outputcellsize(39)=20.10 outputcellsize(40)=20.36 outputcellsize(41)=20.61 outputcellsize(42)=20.86 outputcellsize(43)=21.11 outputcellsize(44)=21.35 outputcellsize(45)=21.59 outputcellsize(46)=21.83 outputcellsize(47)=22.07 outputcellsize(48)=22.30 outputcellsize(49)=22.53 outputcellsize(50)=22.76 outputcellsize(51)=22.99 outputcellsize(52)=23.21 outputcellsize(53)=23.43 outputcellsize(54)=23.65 outputcellsize(55)=23.87 outputcellsize(56)=24.09 outputcellsize(57)=24.30 outputcellsize(58)=24.51 outputcellsize(59)=24.72 outputcellsize(60)=24.93 outputcellsize(61)=25.14 outputcellsize(62)=25.35 outputcellsize(63)=25.55 outputcellsize(64)=25.75 outputcellsize(65)=25.95 outputcellsize(66)=26.15 outputcellsize(67)=26.35 100 outputcellsize(68)=26.54 outputcellsize(69)=26.74 outputcellsize(70)=26.93 'process: Lidar to Grids gp.toolbox = "C:/Model_Builder_Lidar/SensitivityAnalysisTool.tbx" 'for i = 15 to 15 'cellsize=outputcellsize(i) 'cellsize2=outputcellsize2(i) 'iter=iteration(i) 'shape=Shapefile(i) 'watershed=Output_Watershed(i) 'rast=LidarRaster(i) 'gp.LidartoRaster shape, cellsize, cellsize2, rast,iter, lidarpointsamy_mdb 'if i=70 then 'msgbox "loop Lidar to Raster"& i &" finished" 'end if 'next for i = 16 to 69 cellsize=outputcellsize(i) cellsize2=outputcellsize2(i) iter=iteration(i) shape=Shapefile(i) watershed=Output_Watershed(i) rast=LidarRaster(i) gp.RastertoWatershedRaster rast, watershed 'if i=70 then 'msgbox "loop Raster to WatershedRaster "& i &" finished" 'end if next for i = 16 to 69 interwat=Intersect_watershedpoly(i) cellsize=outputcellsize(i) cellsize2=outputcellsize2(i) iter=iteration(i) shape=Shapefile(i) watershed=Output_Watershed(i) rast=LidarRaster(i) watpoly=Output_watershedpoly(i) 101 XS3D= output_XS(i) wat= waterchrono(i) gp.WatershedRastertoPolygonArea watershed, watpoly, interwat, wat if i=14 then msgbox "loop WatershedRaster to watershed Polygon"& i &" finished" end if next 'Derived_Geodatabase = "C:/Model_Builder_Lidar/Output_Sensitivity_Simulation/Watersheds/WatershedPoly.md b" 'WatershedPoly_mdb = "C:/Model_Builder_Lidar/Output_Sensitivity_Simulation/Watersheds/WatershedPoly.md b" 'for i = 68 to 70 'watpoly=Output_watershedpoly(i) 'gp.ToGDB watpoly, WatershedPoly_mdb 'if i=70 then 'msgbox "loop to GDB"& i &" finished" 'end if 'next 102 GLOSSARY ALTM: Airborne Laser Terrain Mapping DEM: Digital Elevation Model DGPS: Differential Global Positioning System DLL: Dynamic Linked Library DTM: Digital Terrain Model FEMA: Federal Emergency Management Agency FIRM: Flood Insurance Rate Map FOV: Field Of View GCP: Ground Control Point GIS: Geographic Information System GPS: Global Positioning System HEC RAS: Hydrologic Engineering Center River Analysis System IDW: Inverse Distance Weighting IFOV: Instantaneous Field of View IMU: Inertial Measurement Unit INS: Inertial Navigation System IRS: Inertial Reference System LIDAR: Light Detection and Ranging POS: Position and Orientation System SFHA: Special Flood Hazard Area TIN: Triangulated Irregular Network 103 BIBLIOGRAPHY Abdalati and Krabill, 1999, ?Calculation of ice velocities in the Jakobshavn Isbrae area using airborne laser altimetry?, Remote Sensing Environ, 67, 194-204 Bates and Marshall, 1999, ?Improving climate model representations of snow hydrology?, Environmental Modeling and Software, 14, 327-334 [Abstract] Bedient, 2003, http://doctorflood.rice.edu/nmd/ accessed in April 2004 Brock et al., 1997; ApJ 329 208 City of Brownsville, RUST Lichter/Jameson, Environment & Infrastructure Consulting Ingineers, Scientist and planners, November 1996, revised in 1997, Flood protection Plan Cobby and Mason, 1999, ?Image processing of airborne scanning laser altimetry for improved river flood modeling?, ISPRS J. Photogrammetry & Remote Sensing, 56, 121- 138 Hodgson et al, 2003, "An Evaluation of LIDAR- and IFSAR-derived Digital Elevation Models in Leaf-on Conditions with USGS Level 1 and Level 2 DEMs", Remote Sensing of Environment 84: 295-308 Kennett and Eiken, 1997, ?Airborne measurement of glacier surface elevation by scanning laser altimeter?, Annals of Glaciology 24, pp. 293-296 Kenward et al, 2000, "Effects of Digital Elevation Model Accuracy on Hydrologic Predictions", Remote Sensing of Environment 74: 432-444 Krabill et al., 2000, ?Greenland ice sheet: high-elevation balance and peripheral thinning?, Science, 289, 428-430 Lohani and Loganathan, 1999, ?Dynamic Modeling of Droughts?, Proceedings of ASCE's Water Resources Planning & Management Conference, Tempe, AZ, June 6-9, 1999 Maidment et al, 2001, http://www.ce.utexas.edu/prof/MAIDMENT/giswr2001/visual/4, accessed in April 2004. Martin and Gutelius, 1997, "Commercial Implications of Topographic Terrain Mapping Using Scanning Airborne Laser Radar"; M. Flood and B. Gutelius, Photogrammetric Engineering and Remote Sensing, April 1997, pp 327-366 104 McCowan, A.D. & Collins, N. (1999), ?The use of Mike 21 for full two-dimensional flood impact assessment?, Users Conference DHI, Horsholm, 8p Wehr and Lohr, 1999, ?Airborne Laser Scanning ? an Introduction and Overview?, ISPRS 54, 1999, pp. 62-82 A. Wehr, U.Lohr, 1999, ASPRS Journal of Photogrammetry & Remote Sensing 54 68- 82. 105 VITA Pierre Gueudet was born in Abidjan, Ivory Coast on April 1979, the son of Isabelle Schoeler and Gilles Gueudet. In June 1999 he completed his DEUG in science of materials at the University of Orsay, France. In 2000 he entered the Ecole Centrale, Lille France, one of the top five ?Grandes Ecoles? of Engineering in France. In 2002, he did an internship at the CNES (National Center of Space Studies) in French Guyana. In August 2002, he entered the Graduate School at the University of Texas at Austin for a M.S in Environmental and Water Resources Engineering. Permanent address: 4 rue Beauregard, 78430 Louveciennes France. This dissertation was typed by the author.