xml version 1.0 encoding UTF-8 standalone no
record xmlns http:www.loc.govMARC21slim xmlns:xsi http:www.w3.org2001XMLSchema-instance xsi:schemaLocation http:www.loc.govstandardsmarcxmlschemaMARC21slim.xsd
leader nam Ka
controlfield tag 001 001790120
007 cr mnu|||uuuuu
008 070613s2006 flu sbm 000 0 eng d
datafield ind1 8 ind2 024
subfield code a E14-SFE0001487
Crosby, David Alexander.
The effect of DEM resolution on the computation of hydrologically significant topographic attributes
h [electronic resource] /
by David Alexander Crosby.
[Tampa, Fla] :
b University of South Florida,
ABSTRACT: Terrain attributes computed from Digital Elevation Models (DEMs) are widely used in hydrology and hydrologic modeling. It is important to consider that the values of the attributes can be different depending on the resolution of the DEM from which they are derived. The question arises as to how much exactly the high-resolution DEMs created through LIDAR remote sensing techniques change the values of the terrain attributes when compared to lower resolution DEMs.In this thesis a LIDAR-derived DEM of 20 feet resolution was resampled using a nearest-neighbour algorithm to various coarser resolutions to examine and quantify the effect of DEM resolution upon a series of hydrologically significant terrain attributes including slope, surface curvature, topographic wetness index, stream power index and stream networks. Values for slope and surface curvature are found to be smaller when computed from lower resolution DEMs; values for the topographic wetness index and stream power index are found to increase as DEM cell size increases.The derived stream networks for each resolution were compared in terms of length per stream order, drainage density, bifurcation ratio, and overall accuracy indicating a loss of small detail, but only a modest change in the overall stream network morphometry. This research suggests that it is possible to establish relationships that quantify the effects of DEM resolution upon hydrologically significant terrain attributes, which can then be considered when processing DEMs from various resolutions for the purpose of parameterizing hydrologic models.
Thesis (M.A.)--University of South Florida, 2006.
Includes bibliographical references.
Text (Electronic thesis) in PDF format.
System requirements: World Wide Web browser and PDF reader.
Mode of access: World Wide Web.
Title from PDF of title page.
Document formatted into pages; contains 135 pages.
Adviser: Paul Zandbergen, Ph.D.
t USF Electronic Theses and Dissertations.
The Effect of DEM Resolution on the Computation of Hydrologically Significant Topographic Attributes by David Alexander Crosby A thesis submitted in partial fulfillment of the requirements for the degree of Master of Arts Department of Geography College of Arts and Sciences University of South Florida Major Professor: Paul Zandbergen, Ph.D. Robert Brinkmann, Ph.D. Mark Ross, Ph.D. Date of Approval: April 6, 2006 Keywords: GIS, LIDAR, DEM, terrain analysis, geomo rphometry Copyright 2006, David Alexander Crosby
Dedication This thesis is dedicated to my family
Acknowledgements I would like to acknowledge and sincerely thank al l of the individuals who have helped me during the research and writing of this t hesis. First, I would like to thank Paul Zandbergen for his guidance, ideas, and expertise. Thank you for the countless hours you spent working through techniques and theory with me and helping make the final product as solid as possible. I would also like to thank m y committee members, Robert Brinkmann and Mark Ross for their numerous suggesti ons and critiques. I truly appreciated working will all of you. I would like to sincerely thank my family for thei r encouragement during my graduate study. I would also like to thank Paul Za ndbergen, Martin Bosman, Steven Reader, Jayajit Chakraborty, and Philip Reeder for helping to provide me with much of the background knowledge and expertise needed to ac complish this research. Thanks as well to my early professors at the University of To ronto for stimulating and piquing my interest in geography and GIS.
i Table of Contents List of Tables iii List of Figures iv Abstract vi Chapter One: Introduction 1 1.1 Background 1 1.2 Goal 4 1.3 Objectives 5 1.4 Description of Study Area 6 Chapter Two: Literature Review 9 2.1 Topographic Models 9 2.2 Creating Terrain Models through Interpolation 11 2.3 LIDAR 12 2.4 Topographic Attributes 17 2.4.1 Primary Topographic Attributes 1 8 2.4.2 Secondary Topographic Attributes 19 2.5 Detailed Description of Attributes of Interest 22 2.5.1 Slope 22 2.5.2 Curvature 22 2.5.3 Upslope Contributing Area and Specific Catc hment Area 23 2.5.4 Topographic Wetness Index 24 2.5.5 Stream Power Index 25 2.6 Hydrologic Modeling and Hydrologically Conditi oned DEMs 25 2.7 Flow Routing Algorithms 27 2.8 TOPMODEL 28 2.9 DEM Resolution 31 2.10 Channel Network Characteristics 32 2.11 Extraction of a Stream Network from a DEM 35 2.12 Threshold for the Extraction of Stream Networ ks 36 Chapter Three: Methodology 39 3.1 Methodology Overview 39 3.2 Data Sources 40 3.3 Processing 41 3.3.1 Resampling Original DEM 41
ii 3.3.2 Computation of Terrain Attributes 41 3.4 Generation of Stream Networks 43 3.5 Comparison of Topographic Attributes 47 3.6 Strahler Stream Ordering of Photogrammetricall y Mapped Streams 48 3.7 Calculation of the Bifurcation Ratio and Strea m Length Ratio 49 3.8 Comparison of Stream Location Agreement using the Kappa Index of Agreement 49 3.9 Vector Analysis of Stream Agreement 51 Chapter Four: Results and Discussion 53 4.1 Topographic Attributes 53 4.1.1 Elevation 55 4.1.2 Slope 60 4.1.3 Overall Curvature 65 4.1.4 Plan Curvature 69 4.1.5 Profile Curvature 73 4.1.6 Stream Power Index 77 4.1.7 Wetness Index 82 4.2 Stream Analysis 88 4.2.1 Stream Network Comparison 8 8 4.2.2 Stream Network Statistical Analysis 100 4.2.3 Comparison of Vector Stream Locations 104 4.2.4 Stream Morphology Statistics 107 Chapter Five: Conclusions 112 References Cited 120 Appendices 125 Appendix A: Geoprocessing Scripts 126
iii List of Tables Table 1: Primary Topographic Attributes that can b e computed by Terrain Analysis from DEM Data, from Wilson and Gallant, 20 00 18 Table 2: Secondary Topographic Attributes that can be computed by Terrain Analysis from DEM Data, from Wilson and Gallant, 20 00 20 Table 3: Number of cells to initiate a stream netw ork 46 Table 4: Summary statistics for elevation 57 Table 5: Results of the Kolmogorov-Smirnov Test on Elevation 59 Table 6: Summary statistics for slope 62 Table 7: Results of the Kolmogorov-Smirnov Test on Slope 64 Table 8: Summary statistics for curvature 67 Table 9: Results of the Kolmogorov-Smirnov Test on the Curvature 69 Table 10: Summary statistics for plan curvature 71 Table 11: Results of the Kolmogorov-Smirnov Test o n the Plan Curvature 73 Table 12: Summary statistics for profile curvature 75 Table 13: Results of the Kolmogorov-Smirnov Test o n the Profile Curvature 77 Table 14: Summary statistics for stream power inde x 79 Table 15: Results of the Kolmogorov-Smirnov Test o n the Stream Power Index 82 Table 16: Summary statistics for wetness index 84 Table 17: Results of the Kolmogorov-Smirnov Test o n the Wetness Index 86 Table 18: Cell Agreement Count Matrix for Boolean Stream Network 101 Table 19: Cell Agreement Count Matrix for 20' DEM Derived vs. Mapped Stream Network 102 Table 20: Kappa Index of Agreement for 20 foot DEM 104 Table 21: Results of Intersection of Mapped Stream Buffer and Derived Stream Networks 105 Table 22: Count of Stream Segments by Order 109 Table 23: Bifurcation Ratios 110 Table 24: Length of Streams by Order (in feet) 110 Table 25: Stream Length Ratios 111
iv List of Figures Figure 1: Study Watershed Showing Topographic Reli ef 7 Figure 2: Study Site in North Carolina 8 Figure 3: LIDAR Data Acquisition, from NC Floodmap s Fact Sheet, 2003 15 Figure 4: Fourth order polynomial in curvature cal culation, from ArcGIS Documentation (ESRI, 2005) 23 Figure 5: Strahler Stream Ordering 34 Figure 6: Image shows road bank impeding derived s treamflow 44 Figure 7: Profile Transect Location in Watershed 5 4 Figure 8: Elevation Variability along a Transect a t Different Cell Sizes 56 Figure 9: Minimum and Maximum Computed Elevation b y Cell Size 58 Figure 10: CFD for Elevation 59 Figure 11: Slope Variability along a Transect at D ifferent Cell Sizes 61 Figure 12: Maximum and Mean Computed Slope by Cell Size 62 Figure 13: Mean Computed Slope by Cell Size with L og Trend Line 63 Figure 14: CFD for Slope 64 Figure 15: Curvature Variability along a Transect at Different Cell Sizes 66 Figure 16: Maximum and Minimum Computed Curvature by Cell Size 67 Figure 17: CFD for Curvature 68 Figure 18: Plan Curvature Variability along a Tran sect at Different Cell Sizes 70 Figure 19: Maximum and Minimum Computed Plan Curva ture by Cell Size 71 Figure 20: CFD for Plan Curvature 72 Figure 21: Profile Curvature Variability along a T ransect at Different Cell Sizes 74 Figure 22: Maximum and Minimum Computed Profile Cu rvature by Cell Size 75 Figure 23: CFD for Profile Curvature 76 Figure 24: Stream Power Index Variability along a Transect at Different Cell Sizes 78 Figure 25: Mean Computed Stream Power Index by Cel l Size 80 Figure 26: CFD for Plan Stream Power Index 81 Figure 27: Wetness Index Variability along a Trans ect at Different Cell Sizes 83 Figure 28: Maximum, Minimum, and Mean Computed Wet ness Index by Cell Size 84 Figure 29: Mean Computed Wetness Index by Cell Siz e with Log Trend Line 85 Figure 30: CFD for Wetness Index 86 Figure 31: Orthophoto Derived Stream Network 88 Figure 32: Streams Derived from 20 foot DEM 89 Figure 33: Streams Derived from 40 foot DEM 90
v Figure 34: Streams Derived from 80 foot DEM 91 Figure 35: Streams Derived from 160 foot DEM 92 Figure 36: Streams Derived from 320 foot DEM 93 Figure 37: Streams Derived from a 20 foot Cell Siz e and Photogrammetrically Mapped Streams 94 Figure 38: Streams Derived from a 40 foot Cell Siz e and Photogrammetrically Mapped Streams 95 Figure 39: Streams Derived from a 80 foot Cell Siz e and Photogrammetrically Mapped Streams 96 Figure 40: Streams Derived from a 160 foot Cell Si ze and Photogrammetrically Mapped Streams 97 Figure 41: Streams Derived from a 320 foot Cell Si ze and Photogrammetrically Mapped Streams 99 Figure 42: All Streams Derived from a DEM and Phot ogrammetrically Mapped Streams 100 Figure 43: Result of Error in Initiating a First-O rder Stream, Cascade Effect 103
vi The Effect of DEM Resolution on the Computation of Hydrologically Significant Topographic Attributes David Alexander Crosby ABSTRACT Terrain attributes computed from Digital Elevation Models (DEMs) are widely used in hydrology and hydrologic modeling. It is i mportant to consider that the values of the attributes can be different depending on the re solution of the DEM from which they are derived. The question arises as to how much ex actly the high-resolution DEMs created through LIDAR remote sensing techniques cha nge the values of the terrain attributes when compared to lower resolution DEMs. In this thesis a LIDAR-derived DEM of 20 feet resol ution was resampled using a nearest-neighbour algorithm to various coarser reso lutions to examine and quantify the effect of DEM resolution upon a series of hydrologi cally significant terrain attributes including slope, surface curvature, topographic wet ness index, stream power index and stream networks. Values for slope and surface curv ature are found to be smaller when computed from lower resolution DEMs; values for the topographic wetness index and stream power index are found to increase as DEM cel l size increases. The derived stream networks for each resolution wer e compared in terms of length per stream order, drainage density, bifurcat ion ratio, and overall accuracy
vii indicating a loss of small detail, but only a modes t change in the overall stream network morphometry. This research suggests that it is pos sible to establish relationships that quantify the effects of DEM resolution upon hydrolo gically significant terrain attributes, which can then be considered when processing DEMs f rom various resolutions for the purpose of parameterizing hydrologic models.
1 Chapter 1 Introduction 1.1 Background Terrain plays a fundamental role in modulating ear th surface and atmospheric processes. Terrain is so fundamental to these proc esses that an understanding of the nature of terrain can also give understanding of th e nature of these processes, in both analytical and subjective terms (Hutchinson and Gal lant, 2000). Knowledge of terrain is important and fundamental to many applications and disciplines, including hydrology, geomorphology, ecology and biology (Wilson and Gall ant, 2000). Terrain can be thought of as the base upon which all surficial ear th processes occur. The following illustrates the linkage between terrain and surfici al earth processes. Surface morphology has an impact both on catchment hydrology and terrain derivatives, such as slope and aspect, and has a di rect impact on solar shading and surface insolation (Wilson and Gallant, 2000). The shape o f the land surface affects the movement and storage of water, nutrients, and other sediment on the landscape. The movement of water and deposition of material in tur n has implications for soil development and geomorphology. Further, this distr ibution of moisture, nutrients, heat, and solar radiation has a direct effect upon photos ynthesizing plants (Wilson and Gallant, 2000). Terrain models generated from elevation dat a, terrain analysis techniques, and
2 Geographic Information Systems (GIS) are tools that allow the informed linkage between surface morphology and biophysical processes. With the advent of new, higher resolution elevation data from new collection metho ds and sensors comes more opportunity to explore there processes at scales no t previously explored. Beven and Moore (1992) noted that digital terrain m odeling is not new in the field of hydrology, and in fact the characterization of t he topography must be important in hydrology. However, what is different now is that the era of GIS, simulation visualization software, and workstation computing p ower has arrived and is rapidly maturing. Further, the quantity of high resolution digital elevation data has expanded dramatically and the revolutionary aspects of the w ork described in their edited volume result from a revolution in the information availab le, rather than from radical changes in methodology. The quantity and quality of digital elevation data has continued to expand, and with the emergence of fine spatial resolution data from sources such as Light Detection and Ranging (LIDAR) sensors, new opportunities for physical-based modeling have been opened up (Atkinson, 2002). Though data quality ma y be higher, it is difficult to determine a priori what level of accuracy is necessary for a given ap plication. LIDAR data represents elevation at a fine spatial resolut ion, but it remains to be seen whether this finer resolution is able to better represent terrai n attributes, and lend better information to the study of hydrology and the predictions generate d by hydrologic models. The importance of elevation data to hydrology has b een clearly articulated in the literature. For example:
3 "Examination of hydrologic processes also requires other information about the surface and subsurface of the earth. Amo ng these, perhaps the most important is topography. Elevation and relate d parameters (slope, aspect, and drainage area) exert an important contr ol on surface and subsurface hydrology and ecosystems. Topography in fluences intercepted radiation, precipitation and runoff movement of sed iment, evaporation, soil moisture, and vegetation characteristics" (Nat ional Research Council, 1991). Terrain analysis work is most often conducted from within a GIS. GIS is sometimes defined as an arrangement of computer har dware, software, and geographic data that people interact with to integrate, analyz e, and visualize the data; identify relationships, patterns, and trends; and find solut ions to problems. Aronoff (1989) defines it simply as a computer-based system used t o store and manipulate geographic information. Longley, et al. (2001) provides an ex cellent summary of the origins of GIS, as well as an explanation of its applicability to d ifferent disciplines. As technology has changed over the years, the range of geoprocessing and analytical tools available within a GIS has grown exponentially (Goodchild, 2004). Lon gley, et al. (2001) has extended the definition to explain that a GIS is a system capabl e of performing any conceivable operation on data obtained from maps. A survey of definitions in the literature, in fact, reveals that data are always mentioned as a central part of the definition. GIS has transformed the way spatial data is analyzed, store d, and manipulated, and encouraged the proliferation of research into data manipulatio n and processing.
4 The relationship between data and GIS tools is mutu ally beneficial. It is the increased availability of high-resolution, continuo us digital elevation data as well as new computerized terrain analysis tools that increased the popularity of research at what Wilson and Gallant (2000) term the toposcale It is the influence of surface morphology at the toposcale that controls catchment hydrology, as well as slope, aspect, and other terrain derivatives. The application that this thesis is concerned with are the hydrologically significant topographic attributes obtained from di gital elevation data, including the accuracy and density of stream networks derived fro m digital elevation data. Scale in the context of GIS and terrain analysis can be though o f as a window of perception or a measuring tool through which a landscape can be vie wed or perceived (Levin cited in Marceau and Hay, 1999). Changing the scale can cha nge the patterns of reality, which has implications for understanding spatial relation ships, and modeling natural systems including a terrain surface (Marceau and Hay, 1999) The specific point of inquiry relevant to this research concerns the effect of DE M scale and grid resolution upon hydrologically related terrain derivatives and stre am networks. What can we gain, if anything from the proliferation of high-resolution digital elevation data acquisition that is poised to take place with the advent of new data co llection methods? 1.2 Goal To determine the effects of DEM resolution upon hyd rologically significant terrain attributes, including stream networks.
5 1.3 Objectives Objective One: To determine the effects of DEM resolution on hyd rologically significant terrain derivatives including slope, pl an curvature, profile curvature, overall curvature, topographic wetness index, and stream po wer index. Hypothesis: Slope: As DEM cell size increases, derived slope values wi ll decrease, and the mean slope in a watershed will decrease. Plan, profile, and overall curvature: As DEM cell size increases, the range of values for the curvature parameters will b ecome limited, and mean values of the curvature parameters in the wate rshed will decrease. Topographic Wetness Index: As DEM cell size increases, watersheds will be modeled to be Â“wetterÂ”, and the mean of the topo graphic index will increase as cell size does. Higher values of the w etness index will be computed, and low values will be less likely. Stream Power Index: As DEM cell size increases, mean stream power index values will increase, as specific catchment a reas will be larger, and the stream power index will increase proportionally to them. These terrain attributes are chosen as they are imp ortant to the study of hydrology and to hydrologic modeling. These attributes of a surface have the ability to influence surficial processes such as flow acceleration and deceleratio n, flow convergence and divergence, soil erosion and deposition, and soil saturation. Objective Two: To determine how the resolution of a DEM affects the morphology and accuracy of the streams derived from it. Hypothesis: Stream Length As DEM cell size increases, the drainage density of streams derived from the DEM will decrease, as stre am lengths will decrease when derived from increasingly coarse DEMs due to lack of detail. Bifurcation Ratio. As DEM cell size increases, the bifurcation ratio between one order of stream and the next will remai n constant with no clear increase or decrease, although the values wil l become less reliable.
6 Stream Accuracy. As DEM cell size increases, the accuracy of stream locations derived from it will decrease, as compare d to streams mapped through orthophotography. The ability to understand the hydrologic properties of a watershed is a central component of understanding and predicting hydrologic response within a watershed. Stream morphology affects the watershed response of many h ydrologic processes including precipitation events, and the erosion, deposition, and movement of sediment. These particular stream morphological attributes were cho sen because of their importance to the study of hydrology and hydrologic modeling. 1.4 Description of Study Area The study area is a small watershed in the southea st corner of Wake County North Carolina as shown in Figures 1 and 2. Located in t he Piedmont area of the state, the climate is extremely variable due to the influence of the topographic variety to the west, and the Gulf Stream on the coast. The area is domi nated by sandy loam well drained soils. Relief in the area is moderate, with elevat ion ranging between 173 feet and 518 feet above sea level, with slopes generally between 2 and 6 percent. Hydrologic characteristics include lakes of various sizes, and a well defined dendridic stream network. The area is generally an erosional enviro nment, and the watershed is approximately 68 square miles in size. This area was chosen for a number of reasons. The North Carolina Floodplain Mapping program collects and compiles high-quality LIDAR DEM data for many of the areas of the state. Freely available LIDAR data fo r this area was available already in DEM format. This particular area has additional be nefits including a natural topography
7 of moderate relief making any effects of DEM resolu tion observable, and the size of the area is large enough to produce statistically signi ficant results and have streams of a sufficient number of orders to make meaningful comp arisons. A very large scale photogrammetrically mapped stream network is availa ble for this area from Wake County government. This stream network will enable meaningful comparisons between DEM derived stream networks and photogrammetrically mapped ones. 013,000 6,500 Feet Figure 1: Study Watershed Showing Topographic Reli ef
8 Figure 2: Study Site in North Carolina Study area in North Carolina
9 Chapter 2 Literature Review 2.1 Topographic Models Moore et al. (1991) describe a Digital Elevation Mo del (DEM) as an ordered array of numbers that represent the spatial distrib ution of elevations above some arbitrary datum in the landscape. Traditionally, this data w ould be captured through photogrammetric techniques. Softcopy or hardcopy s tereoplotters would be used with either aerial photographs or satellite imagery to i nterpret and capture elevation information (Carter, 1988). As this type of equipm ent is not readily available to many potential consumers of DEMs, elevation information is often digitized or captured directly through electronic means, from topographic maps displaying elevation contour information. Maps of this nature are created and a re available from the United States Geological Survey (USGS) and others (Moore, et al., 1991). Irrespective of how the source elevation data is captured when creating a D EM, there are a number of different DEM surface representation structures. Each has di fferent qualities, and the data structure chosen will depend on the end application of the DEM. A structure chosen for computing water storage volumes will differ from on e chosen for computing the topographic attributes of a landscape (Moore, et al ., 1991, Schneider, 2001).
10 It is generally accepted that there are three diff erent ways of structuring digital elevation data. These include regular square-grid raster representations, triangulated irregular networks (TINs), and contours of equal el evation (Wilson and Gallant, 2000). Square-grid DEMs are the most widely used data str ucture. They have several advantages, including their simplicity (single elev ation value stored for each independent grid cell), ease of computer implementation, and co mputational efficiency (Moore et al., 1991, Aronoff 1989, Wilson and Gallant 2000). Ther e are disadvantages as well. First, although there are compression methods available, t he size of the grid mesh will often affect the storage requirements, computational effi ciency, and the quality of the results (Moore et al., 1991). Further, square grids cannot easily handle abrupt changes in topography and computed upslope flow paths will ten d to zigzag across the landscape, making them unrealistic (Moore et al., 1991). Triangulated irregular networks (TINs) are the sec ond data structure. A TIN is composed of nodes, edges, triangles (facets), and h ull polygons. These facets consist of planes joining the three adjacent points in the net work and are constructed in such a way that they meet the Delauney criterion, which ensure s that no vertex lies within the interior of a circle constructed around the triangles in the network (Wilson and Gallant, 2000). TINs sample surface-specific points, such as peaks, ridges, and breaks in slope, and form an irregular network of points stored as a set of x, y, and z values together with pointers to their neighbours in the net (Moore et a l., 1991). TINs are advantageous in that they can easily incorporate discontinuities an d may constitute efficient data structures because the density of the triangles can be varied to match the roughness of the terrain, thus making them flexible and efficient (M oore et al., 1991). Grid-based DEMs
11 do not have this flexibility ability to locally adj ust the grid size to the dimensions of topographic land surface features (Garbrecht and Ma rtz, 2000). Contour-based networks consist of isolines represen ting lines of equal elevation on the surface. These lines are stored as x, y coo rdinate pairs of points located along these lines (Moore, et al., 1991). Contours allow a landscape to be readily visualized, as the isolines are closer together in areas of steep topography, and further apart in areas with less relief. Topographic attributes such as slope, aspect, plan curvature, profile curvature, and the specific catchment area can be derived from gri d-based DEMs, TINs, and contours. However, this is done much more efficiently when us ing the grid-based data structure (Moore, et al., 1991). Details regarding terrain a ttributes and their determination will be discussed later. The rest of this thesis will focu s on the more popular and suitable grid based Digital Elevation Model, hereinafter referred to simply as Â‘DEMÂ’. 2.2 Creating Terrain Models through Interpolation A DEM is often interpolated from point and/or line sample data which has a known or estimated elevation value (Kienzle, 2004). There are many mathematical functions that can be used to fit a smooth surface through a set of sample of points with elevations. More popular ones include local and gl obal polynomial trend surfaces, inverse-distance weighted (IDW) interpolation, spli ning, and kriging. All of these are implemented in popular commercial GIS software appl ications such as ArcGIS (ESRI, 2005) and selectively in standalone analysis packag es such as the Terrain Analysis System (TAS) (Lindsay, 2005). Depending on the imp lementation, each one of these
12 interpolation functions have a set of parameters th at need to be specified, and the results of the interpolation can vary significantly from on e another (Kienzle, 2004). Analysis of the root mean square (RMS) error, an indicator of t he accuracy of the interpolated surface, can be large and will differ between inter polation functions (Kienzle, 2004). One algorithm not particularly suitable for DEM in terpolation from irregularly spaced elevation points (i.e. most sample elevation data sets) is the IDW method. This method tends to create unrealistically shaped terra in features, often referred to as Â“BullÂ’s EyesÂ” (Kienzle, 2004). The interpolation algorithm s that offer the best results include the regular spline with tension and ANUDEM (Huchinson, 1989). ANUDEM is based on a thin-plate spline technique, but has added features that contribute to its creation of DEMs with accurate representation of terrain, and suitab ility for hydrologic modeling. 2.3 LIDAR LIDAR is an acronym for Light Detection and Rangin g. This technology can be used to determine distance, determine chemical comp osition of the atmosphere, or for determining the velocity of a moving object. LIDAR operates in much the same way as traditional RADAR, but instead of emitting radio en ergy and awaiting a reflective response, the LIDAR sensor emits short bursts of la ser light towards an object or a surface. As the emitted light interacts with objec ts and surfaces during its travel, it is reflected back towards the sensor. The time it tak es for the light to be reflected back to the sensor reveals the distance of the object or su rface from the sensor. The intensity of the return is also measured. The distance from the light emitting sensor to the object being detected is determined from mathematical equa tions using the speed of light.
13 LIDAR sensors can be ground based (known as Terrest rial Laser Scanning or TLS), airborne (known as Airborne Laser Scanning or ALS), or spaceborne. Airborne LIDAR sensors emit between 5,000 and 50,0 00 laser pulses per second in a scanning array. Figure 3 shows the most popul ar scanner configuration, which moves back and forth sideways relative to the point s measured on the ground. The average point spacing in the cross-flight direction is determined by the scan angle and flying altitude. Further, both flying height and a irspeed determine the average point spacing in the in-flight direction. Conceptually, each laser pulse can be thought of as a cylinder of light which has a diameter and length, as each laser pulse has a number of properties, including the pulse width which is typi cally between 0.5 and 1 meter in diameter, and a pulse length. The pulse length ref ers to the length of time lapse between the time the laser pulse was turned on and off. Th e three basic types of LIDAR sensing include Differential Absorption LIDAR (DIAL), Doppl er LIDAR, and range finding. DIAL LIDAR is used to measure the chemical composi tion and concentration of molecules in the atmosphere. In order to do this t he LIDAR sensor must emit two different wavelengths of light so that one waveleng th can be partially absorbed by the chemical in the atmosphere and one wavelength can b e reflected back to the sensor. Once the light is reflected back, the chemical conc entration can be determined by measuring the difference in the properties of the t wo light wavelengths (Kavaya, 1999). Doppler LIDAR is used to determine the velocity of a moving object or target. Targets can be either hard or an atmospheric target as small as dust particles. By measuring the velocity of dust particles, scientist s are essentially measuring atmospheric wind speed. This information can be used for futur e weather predictions as well as
14 studying current extreme weather or weather anomali es. Doppler LIDAR works by measuring the Doppler shift of the object or dust p articles. When the light is reflected back to the sensor, the wavelength is either shorte r or longer depending on whether the object or dust particle is moving towards or away f rom the sensor, respectively (Kavaya, 1999). Range finding LIDAR is the most basic kind of LIDA R and is the type used to determine terrain elevations. It operates as descr ibed earlier, by emitting laser light, and measuring the return time and intensity of the puls e. There are 3 components to collecting range finding LIDAR data to produce a DE M. First, one needs the LIDAR sensor itself. This will transmit the laser pulses and record the distance information from the reflected light. Second, the aircraft needs to have a differential Global Positioning System (GPS) on board. The GPS on board records th e location of the sensor in 3 dimensions, including the altitude (z), latitude (y ), longitude (x). Locations are recorded continuously at every point while scanning with the LIDAR sensor relative to a GPS base station, or can be corrected post-mission based on observations of a GPS base station located near the study area (NC Floodmaps, 2003). The third component is referred to as the Inertial Measuring Unit (IMU). The IMU is resp onsible for establishing and recording the pitch, yaw, and roll of the sensor, e ffectively measuring the orientation of the sensor about the three dimensions noted earlier (NC Floodmaps, 2003). In effect, the GPS records the location of the aircraft in three d imensions, and the IMU records in which direction the sensor is oriented. When all o f this information is processed, an accurate location can be specified for every return received by the LIDAR sensor.
15 Figure 3: LIDAR Data Acquisition, from NC Floodmap s Fact Sheet, 2003 LIDAR sensors are capable of receiving multiple re turns, some up to five returns per pulse. This means that a 30-KHz sensor (30,000 pulses per second) must be capable of recording up to 150,000 returns per second. The "first return" recorded by a LIDAR sensor is the first object hit by a laser pulse. T his could be a treetop, rooftop, ground point, or even a suspended object such as a bird in flight or a power transmission line. When a laser pulse hits a soft target (i.e., a fore st canopy), the first return represents the top of that feature. However, a portion of the las er light beam might continue downward below the soft target and hit a tree branch. This would provide a second return, although multiple returns are possible. Theoretically, the last return represents the bare earth
16 surface, but this is sometimes not the case. On oc casion, and in some locations, vegetation is so thick that no portion of the laser pulse can penetrate to the ground. This could be the case with sawgrass, mangrove, and dens e forests where a person on the ground cannot see the sky through the canopy (NC Fl oodmaps, 2003). As LIDAR sensors record multiple returns from each beam of light emitted from the sensor, these multiple returns can be useful fo r differentiating and characterizing terrain, buildings, and even tree canopies. While the first and the last are often the returns of greatest interest, there is ongoing work to determine the utility of the intermediate returns occurring between the first an d last returns. The accuracy of each LIDAR data point depends on m any factors including the flight altitude, the precision of the positioning i nstruments (the GPS and IMU), and the amount and quality of ground control applied. The resolution of the dataset and accuracy of applications built using the dataset depend on t his accuracy as well as on the density of points within the surveyed area. In general terms, the accuracy and density of individual return points increase as the altitude of the fligh t and its speed is decreased. Increasing or decreasing the pulse rate of the laser can also aff ect the point density, or post spacing (Plaster, 2002). The quality of the information re sulting from LIDAR data collection is also dependent on the algorithms that are used to f ilter out vegetation and other ancillary information. Care must be taken when filtering out the data for the last return (terrain) bare earth surface. Each LIDAR data vendor has a p roprietary algorithm to do this, and refining and optimizing this methodology is an ongo ing area of research.
17 2.4 Topographic Attributes Most of the common topographic attributes such as s lope, aspect, plan and profile curvature, and specific catchment area can be deriv ed from elevation data represented as a DEM, TIN, or contours for each and every element as a function of its surroundings (Moore et al., 1991). Individual terrain analysis tools have been classified in various ways based on the characteristics of the computed a ttributes and/or their spatial extent. Some distinguish tools that perform operations on Â“ localÂ” neighbourhoods from those that perform operations on extended or global neigh bourhoods (calculation of upslope drainage areas, etc.) (Burrough and McDonnell, 1998 ). Primary attributes that are computed directly from the DEM are usually distingu ished from secondary or compound attributes that involve combinations of primary att ributes. These usually constitute physically based or empirically derived indices tha t can characterize the spatial variability of specific processes occurring in the landscape (Moore et al., 1991). This same logic is adopted here. Primary attributes inc lude slope, aspect, plan and profile curvature, flow-path length, and upslope contributi ng area. Most of these topographic attributes are calculated from the directional deri vatives of a topographic surface. They can be computed directly with a second-order finite difference scheme or by fitting a bivariate interpolation function z = f(x, y) to the DEM and then calculating the derivatives of the function (Wilson and Gallant, 2000). In man y cases, it may be preferable to calculate a depressionless (hydrologically correcte d) DEM first; specifying one or more rules to determine drainage directions and the conn ectivity of individual elements in order to calculate flow path lengths and upslope co ntributing area. The overall aim is to be able to use the computed attributes to describe the structure, form, catchment position,
18 and surface attributes of hillslopes and stream cha nnels comprising drainage basins. Many researchers, including Band (1986), and Jenson and Domingue (1988), have used computed topographic attributes to generate formal landform classifications (Wilson and Gallant, 2000). The secondary or compound attributes that are compu ted from two or more primary attributes are important because they offer an opportunity to describe pattern as a function of process. Those attributes that quantif y the role played by topography in redistributing water in the landscape and in modify ing the amount of solar radiation received at the surface have important hydrological geomorphological, and ecological consequences in many landscapes. These attributes may affect soil characteristics, as soil is affected by the way water moves through the envi ronment in many landscapes, the distribution and abundance of soil water, the susce ptibility of landscapes to erosion by water, and the distribution and abundance of flora and fauna (Wilson and Gallant, 2000). 2.4.1 Primary Topographic Attributes Table 1: Primary Topographic Attributes that can b e computed by Terrain Analysis from DEM Data, from Wilson and Gallant, 2000 Attribute Definition Significance Altitude Elevation Climate, vegetation, potential energy Upslope height Mean height of upslope area Potential energy Aspect Slope azimuth Solar insolation, evapotranspiration, flora and fau na distribution and abundance Slope Gradient Overland and subsurface flow velocity and runoff rate, precipitation, vegetation, geomorphology, soi l water content, land capability class Upslope slope Mean slope of upslope area Runoff velocity Dispersal Mean slope of dispersal Rate of soil drai nage
19 slope area Catchment slope Average slope over the catchment Time of concentration Upslope area Catchment area above a short length of contour Runoff volume, steady-state runoff rate Catchment area Area draining to catchment outlet Runoff volume Specific catchment area Upslope area per unit width of contour Runoff volume, steady-state runoff rate, soil characteristics, soil water content, geomorphology Flow path length Maximum distance of water flow to a point in the catchment Erosion rates, sediment yield, time of concentratio n Upslope length Mean length of flow paths to a point in the catchment Flow acceleration, erosion rates Dispersal length Distance from a point in the catchment to the outlet Impedance of soil drainage Catchment length Distance from highest point to outlet Overland flow attenuation Profile curvature Slope profile curvature Flow acceleration, erosion/deposition rate, geomorphology Plan curvature Contour curvature Converging/diverging flow, soil water content, soil characteristics Tangential curvature Plan curvature multiplied by slope Provides alternative measure of local flow convergence and divergence Elevation percentile Proportion of cells in a user-defined circle lower than the center cell Relative landscape position, flora and fauna distribution and abundance 2.4.2 Secondary Topographic Attributes A number of the terms which are used in the second ary topographic attribute formulas warrant explanation. A s is the specific catchment area in units of (m2 m-1). The specific catchment area is defined as the surface a rea contributing flow across a width of contour. T is the soil transmissivity when the soil profile is saturated. is the local slope gradient, measured in degrees. A e is the effective specific catchment area. A e differs from A s in that instead of representing all of the area tha t could be expected to contribute
20 to runoff as A s does, A e removes areas of natural storage such as lakes or m arsh areas that would only contribute to runoff during a more significant (1 or 2 year) precipitation event. Table 2: Secondary Topographic Attributes that can be computed by Terrain Analysis from DEM Data, from Wilson and Gallant, 2000 Attribute Definition Significance Topographic wetness index WT = ln tan s T A This equation assumes steady state conditions and describes the spatial distribution and extent of zones of saturation (i.e., variable source areas) for runoff generation as a function of upslope contributing area, soil transmissivity, and slope gradient. W = ln tan s A This particular equation assumes steady state conditions and uniform soil properties (i.e., transmissivity is constant throughout the catchment and equal to unity). This pair of equations predicts zones of saturation where As is large (typically in converging segments of landscapes), is small (at base of concave slopes where slope gradient is reduced), and Ti is small (on shallow soils). These conditions are usually encountered along drainage paths and in zones of water concentration in landscapes. W = ln tan e A This quasi-dynamic index substitutes effective drainage area for upslope contributing area and thereby overcomes
21 limitations of steady-state assumption used in the first pair of equations. Stream-power indices SPI = AS tan Measure of erosive power of flowing water based on assumption that discharge ( q ) is proportional to specific catchment area ( As). Predicts net erosion in areas of profile convexity and tangential concavity (flow acceleration and convergence zones) and net deposition in areas of profile concavity (zones of decreasing flow velocity). LS = (0.4+1) 13. 22 As0.4 0896 .0 sinb1.3 This sediment transport capacity index was derived from unit stream power theory and is equivalent to the length-slope factor in the Revised Universal Soil Loss Equation in certain circumstances. Another form of this equation is sometimes used to predict locations of net erosion and net deposition areas. CIT = AS (tan )2 Variation of stream-power index sometimes used to predict the locations of headwaters of first-order streams (i.e., channel initiation).
22 2.5 Detailed Descriptions of Attributes of Interest 2.5.1 Slope Slope measures the rate of change in elevation and is typically measured as a percent rise or in units of degrees. There are a n umber of different methods that are used to compute slope from DEMs. Methods of calculating slope is an area which has received much attention from researchers, and for w hich there are many references in the literature (Jones, 1998, Horn, 1981). Two common m ethods are the finite differences method, and a method for computing slope based on t he maximum downward slope. Conceptually, the finite differences method fits a plane to the elevation values of a 3x3 cell neighbourhood around the centre cell. The slo pe for the cell is computed from the 3x3 neighbourhood using an average maximum method ( Burrough, 1986). The maximum downward slope method calculates the slope value based on the slope between the centre processing cell, and its lowest neighbou r (Jones, 1998). Slope is significant in hydrology and geomorphology, since it has such an i nfluence on the flux of sediment and water in a landscape (Wilson and Gallant, 2000). 2.5.2 Curvature The curvature of a surface refers to the second der ivative of a surface, or essentially the Â“slope of the slopeÂ”. To calculate the curvature, a fourth order polynomial is fit through the centre processing cell of a 3 x 3 window as shown in Figure 4.
23 Figure 4: Fourth order polynomial in curvature cal culation, from ArcGIS Documentation (ESRI, 2005) Of particular interest are the plan, profile, and o verall curvature values. The plan curvature is the rate of change in aspect, or the r ate of change across a slope. This measure is important for investigating flow converg ence and divergence on a hillslope, and is of great importance in hydrology. The profi le curvature measures the rate of change in slope in the downslope direction. This a ttribute is important for the rate of acceleration of flow, and studies of erosion and de position potential on a surface, and is thus of interest in hydrology and geomorphology. T otal or overall curvature measures the total curvature of a surface. Where the plan a nd profile curvature values measure the curvature of the surface in either an across slope or downslope direction, direction is not a factor in total curvature. Total curvature is also an important attribute for the modeling of flow characteristics (Wilson and Gallant, 2000). 2.5.3 Upslope Contributing Area and Specific Catchm ent Area The upslope contributing area and specific catchmen t area are both important hydrologically related topographic attributes. The upslope contributing area is the area
24 upslope of a length of contour or pixel that would contribute to flow through or over that length of contour or pixel. Closely related to the upslope contributing area is the specific catchment area (SCA). The SCA is different in that it is ratio of the contributing area to a length of contour. In the case of a pixel, the SCA is the ratio of the contributing area to the width of a pixel (Wilson and Gallant, 2000). I n DEM terms, the contributing area can simply be though of the total area of cells which d rain to a given cell. 2.5.4 Topographic Wetness Index The topographic wetness index (Beven and Kirkby, 1 979) is a fundamental component of some hydrologic models, as it measures the potential for soil saturation. It is calculated as the natural logarithm of the speci fic catchment area divided by the tangent of the local slope in degrees (Wilson and G allant, 2000). Topographic Wetness Index = ln tan s A Where A s is the specific catchment area, and is the local slope in degrees This attribute does have a number of assumptions wi th it, including (from Wilson and Gallant, 2000): The steady-state downslope subsurface discharge is the product of average recharge and the specific catchment area The local hydraulic gradient can be approximated by the local slope The saturated hydraulic conductivity of the soil is an exponential function of depth Steady-state conditions are assumed
25 Soil properties, specifically transmissivity are sp atially uniform Other locations in a catchment with the same index value will have the same relationship between the local depth to the water t able and mean depth Other areas of a catchment with the same index valu e will respond in a hydrologically similar way, given the same inputs 2.5.5 Stream Power Index The stream power index is a measure of the erosive force of flowing water. Its components include the unit weight of water (which is a constant), the discharge of water per unit width (similar to specific catchment area) and the local slope. SPI = AS tan Where A s is the specific catchment area, and is the local slope in degrees This measure has been used in the past to predict e phemeral gullies, as well as to model where landscape hardening or conservation measures which would reduce soil erosion should be installed. Variations on this index have been used to predict the locations of the headwaters, or channel initiation location, of first-order streams (Wilson and Gallant, 2000). 2.6 Hydrologic Modeling and Hydrologically Conditioned DEMS Generally, hydrologic modeling with a DEM is a mult i-step process. The reasons for hydrologic modeling with a DEM are varied. In some cases, a hydrologist would be interested in the parameterization of an ecohydrolo gical model such as RHESSys (Creed et al., 1998). Others may be interested in rainfal l-runoff or flood forecasting models, or
26 engineering type applications. Depending on the en d use and result desired, a researcher may take different steps. Typically, the first step in conditioning a DEM for hydrologic modeling is the filling of depressions and pits in the DEM, creatin g a depressionless DEM. This is achieved by raising the values of cells in the depr essions, sinks, or pits, to be equal in elevation with the depressionsÂ’ spill point (Jenson 1992). Other methods used to remove depressions include depression breaching and combin ation methods, such as the Impact Reduction Approach which seeks to fill or breach de pressions depending on which has the least impact to the original DEM (Lindsay and C reed, 2005). Depressions in a DEM can be real surface features, such as lakes, ponds, or reservoirs, or errors in the DEM resulting from interpolation or other factors. At this point, the DEM can be conditioned further to Â‘encodeÂ’ the locations of known streams, or to integrate a vector stream network into a DEM. This process is generally refe rred to as stream burning (Saunders and Maidment, 1996). At the simplest, the DEM valu es are lowered or dropped where streams occur; helping to ensure that flow is route d into these cells during flow accumulation calculations. Another method is known as the AGREE method (Hellweger, 1997), where the DEM elevation value is smoothly lo wered rather than abruptly lowered. This algorithm uses three important reconditioning parameters. The first is known as the vector buffer, and the v alue is entered in terms of number of cells. This represents the number of cel ls around the linear feature class for which the smoothing will occur. The second paramet er is the smooth drop/raise. This is the amount, in vertical units, that the stream will be dropped (if the number is positive) or the fence extruded (if the number is negative). Th is value is used to interpolate the DEM
27 into the buffered area (between the boundary of the buffer and the dropped/raised vector feature). The third parameter is known as the shar p drop/raise Â– this is the additional amount, in vertical units, that the stream will be dropped (if the number is positive) or the fence extruded (if the number is negative). This h as the effect of additional burning/fencing on top of the smooth buffer interpo lation. This parameter is required to ensure the preservation of the linear features used for burning/fencing (Hellweger, 1997). The values that are used for the AGREE reconditioni ng parameters depend on the nature of the DEM itself, and the issues for which resolution is being sought. Often, a trial and error approach is needed before satisfact ory results are obtained (Hellweger, 1997). This integration of vector data into DEMs i s an active area of research, and there are many different ways it can be accomplished, oft en having different results. As a typical second step, flow directions are com puted for each cell in the conditioned DEM. The direction in which water will flow out of a cell into an adjacent one is encoded in the DEM according to a flow-routi ng algorithm (Jenson, 1992). This is a necessary step that needs to occur before a typic al third step, which is the calculation of flow accumulation. This step encodes each cell wit h a number which represents the number of cells upstream of it that would contribut e flow, based on the computed flow directions. Using this raster, and a method for de termining the contributing area threshold, a fully connected stream network can be derived (Jenson, 1992). 2.7 Flow Routing Algorithms The choice of flow routing algorithm in a GIS is on e important component that largely determines the way in which the outflow fro m a given cell will be distributed to
28 one or more downslope cells during the flow accumul ation calculation. A flow routing algorithm functions to transfer flow, including wat er, nutrients, and sediment, to lower (downslope) points or areas in a landscape (Desmet and Govers, 1996). Depending on the flow routing algorithm used to determine the ad jacent downslope areas and points, the calculation of the primary and secondary topographi c derivatives, specific catchment area, stream power index, wetness index, and other topographic indices will change. Flow routing across a surface is an active area of research, and numerous flow routing algorithms have been developed. Of these, the five most commonly used algorithms are the D8 algorithm (OÂ’Callaghan and Mark, 1984), D-In finity (Tarboton, 1997), Rho8 (Fairfield and Leymarie, 1991), FD8 (Quinn et al., 1991), and DEMON (Costa-Cabral and Burges, 1994). Each one of these differs in ho w they route flow downslope, and each will produce a different, albeit sometimes sim ilar, stream network. 2.8 TOPMODEL Topographic attributes, including their scale and a ccuracy, are important inputs to hydrologic models. TOPMODEL is one such model. TO PMODEL is based on the topographic wetness index and depends on it for its modeling where not all hydrologic models do, making the topographic wetness index an important variable in hydrologic modeling. TOPMODEL is a variable contributing area conceptual model in which the predominant factors determining the formation of ru noff are represented by the topography of the basin and a negative exponential law linking the transmissivity of the soil with the vertical distance from the ground lev el. In this model the total flow is calculated as the sum of two terms: surface runoff and flow in the saturated zone. Chairat
29 and Delleur (1993) quantified the effects of DEM re solution and contour length on the distribution of the topographic wetness index as us ed by TOPMODEL and the modelÂ’s peak flow predictions. Wolock and Price (1994) and Zhang and Montgomery (1994) also examined the effects of DEM source scale and DEM ce ll spacing on the topographic wetness index and TOPMODEL watershed model predicti ons. TOPMODEL represents catchment topography by means o f a topographic index, ln(As/tan ), where Â‘AsÂ’ is the area draining through a grid square per un it length of contour and Â‘tan Â’ is the average outflow gradient from the square. The index is calculated from a DEM across a grid covering the ca tchment. The grid must be sufficiently fine to resolve important characterist ics and slope formations. A high index value usually indicates a wet part of the catchment ; which can arise either from a large contributing drainage area or from areas with very low slope. Areas with low index values are usually drier, resulting from either ste ep slope or a small contributing drainage area. Grid squares with the same index values are assumed to behave in a hydrologically similar manner. Topography is now recognized as a first-order contr ol on the hydrologic response of a catchment to rainfall. Topography plays an ex tremely important role in determining the catchment scale flow pathways resulting from th e downward force of gravity (Brasington and Richards, 1998). TOPMODEL is a cat chment scale rainfall-runoff model which makes an explicit link between catchmen t topography (and thus DEMs and related topographic attributes) and the generation of streamflow. The model is based on a spatially distributed topog raphic wetness index, a compound topographic attribute that can be computed from DEMs. Apart from
30 TOPMODEL, the topographic wetness index is used ext ensively in hydrology, agriculture, geomorphology and vegetation studies, as it represents the spatial distribution of groundwater recharge areas, surface saturation p otential, soil moisture, and indicates areas with potential for runoff generation (Kienzle 2004). Water will tend to accumulate in areas with a high value of the topographic wetne ss index. DEM resolution is an extremely important factor in this case, since it h as been noted that coarse DEMs tend to model landscapes as being wetter, whereas finer DEM s model landscapes to be drier (Wolock and McCabe, 2000). TOPMODEL results have b een shown to be sensitive to DEM resolution (Brasington and Richards, 1998). P revious studies have shown significant differences in the probability distribu tions of topographic index computed from DEMs of different resolutions; including 12.5m and 50m (Quinn et al., 1991). Zhang and Montgomery (1994) reported that the mean of the topographic index increased progressively with DEM size, for DEM sizes ranging from 4 to 90m. Wolock and Price (1994) made the same conclusion in their study of t he effect that DEM resolution and the scale of the data used to derive the DEM had on the computed topographic index. The topographic wetness index has much influence on TOPMODEL hydrologic modeling results. The effect of DEM resolution has been investigated specifically in terms of TOPMODELÂ’s modeling predictions. Zhang an d Montgomery (1994) found that the shift of the index towards higher values i ncreased the rate of predicted peak streamflow.
31 2.9 DEM Resolution It is well known that the grid cell size of a raste r DEM has a significant effect on the elevation derivatives, such as slope, aspect, c urvature and the topographic wetness index. The values of these derivatives will depend on the DEM resolution (Kienzle, 2004). In one study, Kienzle (2004) determined tha t both slope and soil erosion estimations increase with a decrease in the grid ce ll size used to derive them (Kienzle, 2004). In the study conducted by Saulnier et al. ( 1997), it was found that the topographic wetness index increased as the grid cell size did. Essentially, this means that catchments are modeled to be wetter when using a coarse grid c ell size, and drier when using a smaller grid cell size (Saulnier quoted in Kienzle, 2004). When Zhang and Montgomery (1994) investigated two small catchments in the sta tes of Oregon and California, with grid cell sizes ranging from 2 to 90 m, different v alues for slope were found. Specifically for the two catchments, a mean slope of 65% and 34% was found with the 2 m grid, as compared to 41% and 29% when using the 90 m grid (K ienzle, 2004). Further, as determined by Kienzle (2004), as grid cell size inc reases, the slope values get smaller. DEMÂ’s with grid cell sizes over 25 m are not able t o identify steep slopes successfully. Also, as the grid cell size increases, a considerab le underestimation of slope value occurs. In the study conducted by Brasington and Richards ( 1998), they found that hillslope gradient is sensitive to grid cell size, with progressively fewer cells showing a steep slope, as grid cell size increases. This is indicated through the mean slope. Brasington and Richards (1998) also looked at the e ffects solely on upslope contributing area (represented as a ). This parameter can be defined as the total upsl ope area draining through a single pixel. They found that large grid cell sizes necessarily increase the
32 minimum contributing area (the square of the grid s ize), and increase the values of ln( a ) as grid cell size increases. The impact of increas ing grid cell size in lowering slopes and increasing contributing area is reflected in the co mpound topographic wetness index. As the grid cell size increases, the mean values of th e topographic wetness index increased as well. It is also accepted that different qualities of DEM s can influence a derived stream network, in addition to the stream network algorith m used. Zhou and Liu (2002), quoting Zhang and Montgomery (1994) observe that the spatia l data structure of DEMs such as data precision, grid resolution and orientation aff ect how terrain is analyzed, including the generation of stream networks. One method sugg ested by Zhou and Liu (2001) is the visual comparison of derived topographic informatio n, such as a stream network, against a Â‘common knowledgeÂ’ of what would be expected, suc h as streams derived from a high quality orthophoto. Statistical comparison between observed and modeled attributes is also suggested. DEMs used in hydrologic analysis have to be of suff icient resolution to capture the variability of the terrain that they represent, as it is the terrain that plays such an efficient role in the determination of landscape fl ow pathways (Brasington and Richards, 1998). An increase in inexpensive computing power and hard drive storage has encouraged many to use DEMs of high resolution. 2.10 Channel Network Characteristics An important aspect of this work is the definition of channel network characteristics. Stream systems are a form of a ne twork, having nodes and links and
33 branches that connect all of the parts together. N etworks in general can be characterized and analyzed with respect to two main sets of prope rties. These are the topologic aspects of the network and the geometric aspects of the net work. The topology of the network refers to its interconnectedness. The more relevan t properties in this case are the geometric properties, which can include the length, shape, area, relief, orientation, and arrangement of the streams (Summerfield, 1991). The most basic component of a stream network is an individual stream segment or link. A stream segment is simply a segment of a st ream between two channel junctions, or in the case of a first-order stream, is a segmen t where flow into a network originates (i.e. it has no tributaries) (Summerfield, 1991). This concept of stream order is also of vital importance to the study of basin morphology, description of stream networks, and to applications concerning hydrologic modeling as it c an be related to the relative discharge of an individual channel segment. There have been a number of stream ordering techni ques proposed by hydrologists and geomorphologists. Of these, the t wo most common methods are known as the Strahler ordering method, and the Shreve ord ering method, with the most popular one being the Strahler method (Lindsay, 2005). The Strahler ordering method designates all stream network headwater tributaries as a first -order stream. A second-order segment is formed where two first-order segments are joined joining two second-order segments produce a third-order segment, and so on. As a pro perty of the Strahler ordering method, there is no increase in stream order when a particu lar stream segment is joined with a stream segment with a lower order. The Strahler or dering system has been found to be statistically related to important basin morphologi cal properties in a wide range of
34 environments (Summerfield, 1991). Figure 5 shows t he graphic representation of a Strahler ordered stream network. Figure 5: Strahler Stream Ordering The drainage (stream) density is a geometric prope rty of a stream network, and is simply defined as the mean length of stream channel s per unit area, functionally the length of streams in a watershed, divided by the wa tershed area (Summerfield, 1991). The drainage density is one of the most important m easures of basin morphometry as well as an important indicator of the relation betw een climate, vegetation, and the resistance of rock and soil to erosion. The draina ge density gives an idea of the resistance of the land surface to erosion by flowin g water, and is related to climate and lithology. Climate is a natural factor influencing observed drainage densities. Very high
35 drainage densities can be observed in semi-arid env ironments, which are a consequence of prevalent surface runoff, and the ease with whic h a new channel can be originated. The stream order is an important concept that is u sed to compute other stream morphology characteristics. One of these of import ance to this study is a topographic metric termed the bifurcation ratio. The bifurcati on ratio is defined as the ratio of the number of streams of an order to the number of stre ams of the next highest order, and is calculated as: Rb = Ni/Ni+1. The number generally varies somewhat between diffe rent successive orders in a watershed (i.e. between firs t and second order, and second and third order). For this reason, a mean bifurcation ratio is usually reported for a watershed to enable comparisons. Typical values for the bifu rcation ratio range between 3 and 5, with 3 considered to be theoretically normal. Bifu rcation ratios as high as 10 can be observed in highly elongated watersheds with certai n combinations of lithologic properties. An example would be a watershed with a lternating outcrops of relatively hard and soft lithologies (Summerfield, 1991). All metrics related to stream network topology an d geometric properties are affected by the accuracy of the streams used to der ive them. Stream accuracy is the relationship between the measured or observed strea m network and the actual location of the stream network in the landscape. 2.11 Extraction of a Stream Network from a DEM There is a well developed approach to hydrologic an alysis using a GIS. The foundational component of all analysis is a DEM. T he first step in conducting an analysis is creating a depressionless DEM, ensuring that all flow routed with the GIS is
36 directed to a downstream outlet, and is not Â“stuckÂ” in an artificial (or natural) depression. A number of algorithms can be used to create the de pressionless DEM, which is usually achieved by filling depressions or pits (Doan, 2000 ). Typically, when there is stream or Â“blue lineÂ” data of a high accuracy available, this data will be used in conjunction with the DEM to aid in the extraction of a digital strea m network. There are numerous methods and algorithms used to accomplish this, non e of which will be used in this work. The reason for the exclusion of this step is that t he DEM itself will be evaluated for utility in digital stream network extraction, and i ntroducing a conditioning such as this would obviously bias the result. From the DEM, rasters of flow direction and flow ac cumulation can be derived. The flow direction grid partitions the surface, def ining how flow, sediment, nutrients and other constituents flow over the surface from one r aster cell to adjacent ones. Different algorithms exist for routing this flow, the most co mmon being the Deterministic 8 (D-8) method. Flow accumulation can be determined using the flow direction raster, indicating how many Â“upslopeÂ” pixels drain through a single do wnslope pixel. Much methodological and empirical research has been cond ucted, and is being conducted now, to determine how best to extract stream networks an d topographic attributes from DEMs (Tarboton, 1991). These topographic attributes and techniques form the basis for most GIS-based hydrologic modeling. 2.12 Threshold for the Extraction of Stream Network s Numerous examples of techniques for digital channe l extraction from DEMs can be found in the literature. Peuker and Douglas (19 75) described a method of extracting
37 stream points from DEMs by flagging the pixel with the highest elevation in a window of four cells. After one complete pass of the moving window through the DEM the pixels that remain unflagged are determined to be drainage courses. One problem with technique however, is that the pixels representing the drainage courses are not necessarily connected. Band (1986) describes several improveme nts to this algorithm including procedures to thin and connect these potentially no n-contiguous pixels. Band (1986) also described a method of flagging upwardly concave pix els in a three by three moving window, which is a change the Peuker and Douglas (1 975) algorithm. While both of these methods result in a potential drainage channel network, there is no physical basis for the network. The work of OÂ’Callaghan and Mark (1984) suggests a physical basis for channel extraction. Their dig ital extraction method defines channel networks as all pixels with an accumulated area abo ve a threshold. That is the accumulated area of all Â“upslopeÂ” pixels. This ide a is classified as hydrologic, and works on the assumption that: Â“Drainage represents those points at which runoff i s sufficiently concentrated that fluvial processes dominate over s lope processes. If the spatial concentration of surface runoff is simulate d, then those points at which this runoff exceeds some threshold can be con sidered to be the drainage networkÂ” (Mark quoted in Sole and Valanzan o, 1996). The drainage network is found by determining a drai nage direction matrix, identifying and removing sinks in the surface, defining a weigh t matrix, calculating an accumulated area matrix (based on a flow accumulation and flow routing algorithm) and selecting those pixels with a contributing area above a given threshold (OÂ’Callaghan and Mark,
38 1984). Their work does not give any indication as to what an appropriate, or physically justifiable contributing area is though. Tarboton (1991) proposes methods for extracting di gital channel networks from DEMs based on a physically justifiable accumulation area threshold. While considering that digital channels are most appropriately extrac ted from DEMs at a scale that is appropriate, perhaps the easiest way of determining the cumulative accumulation threshold is to compare the digital network extract ed at different thresholds to the Â“blue linesÂ” drawn on traditional paper maps (Tarboton, 1 991). Further, one could compare channels extracted with different accumulation area thresholds to photogrammetrically mapped streams. There are arbitrary decisions made with both methods though, especially considering that most first-order stream s are ephemeral and present only during the Â“wet seasonÂ” (Strahler, 1957). Tarboton (1991) proposes two methods for extractin g the most dense or highest resolution drainage networks quantitatively, by see king to satisfy scaling laws that have been found to hold for channel networks. These law s include a constant drop property and a power law which scales slope with area.
39 Chapter 3 Methodology 3.1 Methodology Overview This chapter will describe the data sources and met hods employed to answer the two research objectives specified in the first chap ter. In brief, tiled raster DEM data generated from LIDAR points were acquired and merge d for a watershed using ArcGIS desktop GIS software. The resulting DEM had a reso lution of 20 feet. This DEM was resampled to successively coarser resolutions. Fro m these DEMs, the various topographic attributes to be studied were derived. These included the slope, plan curvature, profile curvature, overall surface curva ture, topographic wetness index, and the stream power index. To facilitate the repetitive o perations, to ensure repeatability, to reduce error, and facilitate making changes should they be needed, geoprocessing models were built using the ArcToolbox in ArcGIS. (Appendix A) Stream networks were derived from the DEM, by comp leting a process of hydrologically conditioning the DEM (ensuring the s urface is conditioned to allow for downslope drainage), including only the filling of pits and depressions. Flow direction and flow accumulation rasters were created for each DEM resolution, from which streams were extracted, based on a contributing are a threshold, chosen to create a stream network that most closely matched the photogrammetr ically mapped stream network.
40 Once the computation of the topographic attributes was complete, and stream delineation was complete and repeated for each DEM, the results could were tabulated and compared. 3.2 Data Sources The first methodological task was to acquire elevat ion data, collected through LIDAR remote sensing. Using LIDAR data would allow the comparison between high resolution data and more popular existing data sets such as United States Geological Survey (USGS) 30m DEMs. This data was available fo r many parts of North Carolina; North Carolina has flown many areas of the state un der the North Carolina Flood Maps program. Also beneficial was the availability of a large-sca le, high accuracy stream network for this study area, in a digital spatial format. The stream network data for Wake County, North Carolina met the National Map Accuracy Standa rds (NMAS) for 1:1200 data (+/5 ft). It is current to June 20, 2005, and was compi led by Surdex, Inc. for Wake County at a scale of 1:1200 from aerial photography at 1" = 6 00'. The availability of existing mapped data allowed the comparison of stream networ ks extracted from the DEMs, and those mapped from large scale orthophotography.
41 3.3 Processing 3.3.1 Resampling Original DEM The original LIDAR DEM acquired for the study area was resampled to various coarser resolutions, including 40 feet, 80 feet, 16 0 feet, 320 feet and 640 feet. To effectively resample the DEMs successively to the c oarser resolutions, a nearestneighbour method was employed. The original DEM wa s resampled as follows, following a method used by Wolock and McCabe (2000) : To create a 40 ft DEM, the original DEM was sampled at every second point To create an 80 ft DEM, the original DEM was sample d at every fourth point To create a 160 ft DEM, the original DEM was sample d at every eighth point To create a 320 ft DEM, the original DEM was sample d at every sixteenth point To create a 640 ft DEM, the original DEM was sample d at every thirty second point In all cases, the resampling of the DEM from 20 ft to a coarser resolution resulted in a raster containing a loss of terrain information. 3.3.2 Computation of Terrain Attributes Once the DEM was resampled, the topographic derivat ives were calculated and compared to quantify the change. Calculation of th e first terrain derivative of slope was carried out for each resolution DEM using HornÂ’s me thod (Horn, 1981) employed in ArcGIS software. Degrees of slope were calculated, as op posed to percentages. Values for plan curvature, profile curvature, and curvatur e were calculated in ArcGIS using the
42 CURVATURE function within ArcGIS (Moore et. al., 1991). Each of these terrain attributes were computed for each DEM of various re solution. Compound topographic attributes including the topog raphic wetness index and the stream power index were calculated in the stand alone Terrain Analysis System (TAS) hydrologic and terrain analysis package (Lindsay, 2 005). TAS has the algorithms for calculating these attributes built in, and completi ng this step in TAS rather than ArcGIS reduced the need to manipulate data within ArcGIS. TAS uses its own raster data format, known as a TAS Image File. Generating the TAS Image File was completed by a two step process, including the conversion of the o riginal ArcGIS DEM to an ASCIIformat file and subsequent import of them to TAS us ing the Â“Import RasterÂ” function. The results of the computation were a series of ras ters of the topographic wetness index and stream power index for each of the various rast er resolutions. The TAS image files for wetness index and stream power index were then exported from TAS as an ASCII format using the Â“Export RasterÂ” function, then imp orted to ArcGIS as a raster, using the ASCII to raster function in ArcToolbox. When the six TAS images for the stream power index were converted to ASCII and then to ArcGIS rasters, areas outside the study area acquired val ues of 0. In order to calculate meaningful statistics for each raster, th ese areas with the value of 0 outside of the study area were converted to null values, using the ArcInfo GRID SETNULL command, in the form: gridnameXXnull = setnull(gridnameXX == 0, gridnameX X)
43 This command would evaluate on a cell-by-cell basis whether a cells value was 0. If the expression evaluated to be true, it was given a nul l value, and if the expression evaluated false, it was given the existing value. 3.4 Generation of Stream Networks In order to compare the stream networks that result from DEMs of different resolutions, a steam network was delineated for eac h DEM, and compared to known stream locations derived from large-scale orthophot ography. To create the stream networks, the following method was employed. A wel l documented problem in stream network derivation from LIDAR DEMs is the effect of road banks and anthropogenic modification of the surface upon surface hydrology (Duke et al, 2003). In many cases, a road bank has been built up which serves as a dam, preventing downslope flow. In reality, a culvert or pipe is placed beneath the ro ad bank, allowing surface runoff to maintain a route close to what is natural. However to ensure that a derived flow network matches the natural one as closely as possible, and flow is not intercepted by the road banks, slight DEM modifications were completed to c orrect for this condition. Areas where flow was unnaturally intercepted were i dentified and corrected. First, a hydrologically correct depressionless DEM was computed from the original 20 foot resolution DEM. This was completed by using t he FILL SINKS tool in ArcGIS software. The result was a depressionless DEM, in which there are no artificial or natural pits or sinks, ensuring that all flow in the DEM is in a downslope direction and is directed to an outlet, allowing the creation of a flow direc tion raster. The flow direction raster
44 was computed for this depressionless DEM with the D 8 method (OÂ’Callaghan and Mark, 1984). Using the flow direction raster, a flow accumulatio n raster was computed for the DEM. From the flow accumulation raster a stream ne twork was queried out by using a stream accumulation threshold or critical source th reshold. An arbitrary threshold value of 400 cells to initiate a stream network was selec ted. Using this stream network, the derived flow was visually compared to the photogram metrically mapped stream network, to identify those areas where flow was interrupted by a road bank or other artificial impediment. Figure 6 gives an example of this effe ct. Road centerlines representing road banks are green lines; a derived stream networ k is represented by red lines, and the photogrammetrically mapped stream network is repres ented by blue lines. The blue stream is observed to flow Â“throughÂ” the road bank, presumably diverted through a culvert, where the red lines indicate the derived f low traveling along the road bank, not following where the actual stream location is. Derived Stream Network Mapped Stream Network Roads 04809601,4401,920 240 Feet Figure 6: Image shows road bank impeding derived s treamflow
45 A number of techniques were evaluated for their abi lity to modify the DEM in only those areas identified to have artificially in terrupted flow. Of these, one method was eventually used. One of the attributes of the Wake County streams shapefile (Adler, 2001) was stream type. Of particular interest to t his exercise was the Â‘connectorÂ’ stream type. Connectors in this instance were identified to be artificial stream segments, such as those that are diverted through culverts and under road banks. Attempts were made to Â“burnÂ” these connector streams into the original DE M, lowering the elevation values at these locations, to divert flow to the lower elevat ions. Once the stream burning was completed, a stream network was generated again, th rough the iterative process of filling depressions, deriving flow directions and accumulat ions. The outcome was not as intended though, with the stream network remaining unchanged. Upon careful examination of the DEM and the derived flow network s, it was obvious that a gradient needed to be used to modify the DEM in such a way t hat flow would be directed downslope, essentially through a funnel to a connec tor location. To accomplish this task, the AGREE algorithm (Hellweger, 1997) was used. In order to have the least modification of the DEM possible, a number of diffe rent reconditioning parameters were used in the reconditioning. The AGREE parameters e ventually used to create the reconditioned DEM were as follows: Vector buffer (cells) = 10 Smooth Drop (feet) = 5 Sharp Drop = 5 This step used a total of only 31 connector segment s which allowed the alteration of the DEM only where necessary. Reconditioning the DEM u sing AGREE accomplished the
46 task of only slightly modifying the required areas and leaving the remainder of the DEM unaltered, as well as forcing flow to flow where it should underneath road banks. The 20 foot DEM (DEM20) was reconditioned, resulting in AG REEDEM20. Pits and sinks were filled in this DEM resulting in FILLDEM20. AG REEDEM20 was resampled to 40, 80, 160, 320, and 640 foot resolutions, and each ha d pits and sinks subsequently filled, resulting in FILLDEM40, FILLDEM80, FILLDEM160, FILL DEM320, and FILLDEM640. Each AGREEDEM had sinks and pits fille d, as the resampling process introduces small pits and sinks during its executio n. These DEMs became the ones to be used in the final stream network generation and ana lysis, and watershed delineation. To prevent simply choosing an arbitrary stream accu mulation threshold, a series of stream networks were calculated for FILLDEM20 ba sed on values of a stream accumulation threshold of 200 to 700 cells. The re sulting stream networks were compared visually to determine which one most close ly resembled the photogrammetrically mapped stream network in terms of the initiation point of first order streams. The stream accumulation threshold event ually chosen was 400 cells to initiate a stream network. This equates to 160,000 square f eet, or 1.5 ha. The following table gives the number of cells at each resolution, used to initiate and derive a stream network. Table 3: Number of cells to initiate a stream netw ork Cell Size (ft) Number of cells to initiate a stream Corresponding square feet Hectares 20 400 160000 1.49 40 100 160000 1.49 80 25 160000 1.49 160 6 153600 1.43 320 2 204800 1.90
47 Following the generation of the stream network, the Strahler stream order was assigned to the stream links in the network through the STREAMORDER tool present in ArcGIS software, for each resolution. The resulting stre am network was then converted to the vector stream features in ArcGIS using the Streams to Feature tool (without simplification) so that this data was available for subsequent analysis. 3.5 Comparison of Topographic Attributes A number of descriptive, quantitative and qualitat ive results were compiled and reported. For the topographic attributes slope, pl an curvature, profile curvature, surface curvature, topographic wetness index, and the strea m power index, descriptive statistics including mean, maximum, minimum, and standard devi ation were calculated using ArcGIS and reported in tabular format. The Cumulative Fr equency Distributions (CFDs) of each of these attributes were graphed to observe any visual difference in the distributions. To enable a visual comparison of the effect of res olution on terrain attributes, profile charts were created to show how the values of a terrain attribute vary over a twomile transect of the watershed. First, the locatio n of a two-mile transect was identified that included variation in the terrain observed thr ough inspection of the 20 foot DEM with a hill shaded rendering applied. The location of the transect on the surface is shown in Figure 7. The LIDAR Data Handler ArcMap extension (NOAA) was then used to extract the terrain attribute value every 20 feet f rom every raster for all attributes and resolutions. This data was then exported to a comm ercial spreadsheet software package for creation of the profile charts.
48 To enable a more statistical comparison, the stati stical distributions of each of the topographic attributes were compared using a statis tical test known as the KolmogorovSmirnov (K-S) test. This statistical test is used to determine if there is a significant difference between two CFDs. The test was conducte d on the CFDs for each topographic attribute at different resolutions, to compare and determine if there is a significant difference in the distributions. The tests were co nducted using the Terrain Analysis System (TAS) software. 3.6 Strahler Stream Ordering of Photogrammetrically Mapped Streams While ArcGIS has the facility to attribute a raster stream netw ork with Strahler stream order information when derived from an avail able DEM, the base software product cannot accomplish the ordering of vector st ream networks without customization. In order to compare, by Strahler order, the streams that are delineated from the DEM surface in the GIS, and those vectors of streams ac quired through photogrammetric mapping methods, the vector streams need to have an order value assigned to them. This task was completed using a customized software prod uct named RivEx ( http://www.rivex.co.uk ). Some minor modifications to the stream network were required to allow the software to correctly identif y the stream order. First, all stream segments were checked to ensure that they flow in t he proper downstream direction. The original stream network was also constructed in suc h a way that there could be several physical stream segments that comprise one actual s tream. This could cause the count of stream segments to be wrong. To avoid this problem all pseudo-nodes were removed from the original stream network. This ensured tha t any given physical stream segment
49 was represented by only one line in the GIS. The o rdered streams were then used to compare the spatial location and orders of analytic ally delineated streams. 3.7 Calculation of the Bifurcation Ratio and Stream Length Ratio The calculation of the bifurcation ratio was compl eted for streams delineated from DEMs of the various resolutions as well as the phot ogrammetrically mapped stream network. The number of stream segments in each str eam order for each stream network was summarized in ArcGIS software and the bifurcation ratios were manually calculated. The lengths of all streams by order we re summarized in ArcGIS for each stream network giving a total of all first-order th rough sixth-order streams, for all networks. The stream length ratios were then calcu lated manually (Tarboton, 1996). 3.8 Comparison of Stream Location Agreement using t he Kappa Index of Agreement The Kappa Index of Agreement (KIA) is a statistical test that can be used to test raster independence, and to quantify the similarity (or dissimilarity) of rasters. The KIA is a number between 0 and 1, with a higher number i ndicating more agreement (Chuang, 2001). The KIA was used to compare the stream netw ork generated through automated delineation from the 20 foot DEM, with the photogra mmetrically mapped stream network (Melville and Martz, 2004). To facilitate this, th e rasters being compared need to have the same resolution, and need to have the same orig in and alignment. The original vector photogrammetrically mapped stream network was conve rted to a raster format using a cell size of 20 feet, using the 20 foot DEM as a ba se to snap to.
50 The KIA was determined for overall raster agreemen t and for agreement by order through processing in ArcGIS. Essentially, a cross-tabulation of the rasters w as performed to generate the KIA statistic. Each rast er was converted to a Boolean stream network raster, where the value of 1 represented a stream cell, and 0 represented a nonstream cell. Through a map algebra combine operati on, the stream agreement could be determined. The output of the operation is a matri x showing the four potential combinations of the agreement, including 0-0 (not s treams in both stream rasters), 1-0 or 0-1 (a stream in one raster but not the other) and 1-1 (a stream in both rasters indicating agreement). Using this matrix, the KIA was determi ned. To facilitate the generation of KIA by order, a ma p algebra 'combinatorial and' operation was completed with the photogrammetricall y mapped streams raster and the streams derived from the 20 foot DEM. This operati on resulted in a table showing giving a unique value to cell overlap, by order. This tab le provided the data to generate a matrix showing stream agreement by order, and the subseque nt generation of the KIA for each order of stream. This same matrix was used to gene rate the KIA for the derived and mapped stream network considering all orders. The KIA was only calculated for the results of the stream network generated from the 20 foot DEM. Although the KIA is a flexible st atistic used in a great number of areas and disciplines, in this context the KIA only works on raster data and using a single cell size in the analysis is a requirement. This makes it impractical to compare the results of a stream network generated from cell sizes of many re solutions (i.e. 40, 80, 160, and 320) to the rasterized photogrammetrically mapped stream network, with a cell size of 20 feet.
51 Although the KIA is a very good measure of agreemen t, another method was required to compare the agreement across all resolutions. 3.9 Vector Analysis of Stream Agreement In order to facilitate the comparison of stream net works delineated from rasters with different resolutions, and the streams that we re photogrammetrically mapped, the different networks were compared quantitatively. T here were a number of steps involved in the processing to accomplish this. First, the p hotogrammetrically mapped streams were ordered using the RivEx stream ordering tool, which assigned each unique stream arc a Strahler order attribute. Each of the raster stream networks that were derived from the DEMs of increasingly coarse resolutions were co nverted to vector stream networks through processing in ArcGIS. Then, the photogrammetrically mapped streams wer e buffered by 20 feet, the smallest DEM cell size use d in the study, and the buffers were concurrently dissolved based on the stream order at tribute. There was a small overlap of the buffers that occurred where stream segments con nect, meaning that in a very few areas, it is possible for a stream segment to be lo cated in more than one buffer, and attributed to the buffer of a stream network to whi ch it does not belong. Each of the vector stream networks at the various resolutions w ere then intersected with the buffered photogrammetrically mapped streams through a proces s similar to that employed by Kenny and Matthews (2005). This process assigned e ach stream segment that intersected the buffer with the length within and the order of the buffer it intersected. By determining where the Strahler stream attribute fro m the vectorized stream networks matched the Strahler attribute from the buffered st reams, a matrix of agreement was
52 constructed. This matrix shows the lengths by orde r of the stream segments that were intercepted correctly. This matrix also supports t he calculation of the percentage, by length, of correct interception.
53 Chapter 4 Results and Discussion 4.1 Topographic Attributes For all the terrain attributes, descriptive statis tics were generated for the extent of the study watershed. Depending on the attribute, d ifferent statistics were graphed to visually show the trend. There was no general patt ern found that was consistent for all derivatives. Instead, each computed derivative sho wed a different resolution dependency. In order to visualize how the values for computed terrain derivatives change over the surface at different cell sizes profile graphs were created for each attribute. Figure 7 shows the location of the transect on the surface a long which each attribute was extracted.
54 06,40012,800 3,200Feet Figure 7: Profile Transect Location in Watershed For each topographic attribute a cumulative distri bution function was created. The cumulative frequency distributions (CFDs) show each topographic attributeÂ’s values, at each DEM resolution, as a percentage of the terr ain surface. In general, the cumulative
55 frequency diagrams are excellent indicators of chan ges in the distribution of the various terrain attributes at different resolutions. The d iagrams are further indication of the direction of change (if any) in the values at each cell size. In some cases, there is an expected parallel shift in the CFD graph for a terr ain attribute such as the topographic wetness index. For others, like the curvature para meters where the curvature values are expected to become limited when computed from an in creased cell size, the diagrams will show high representation at extreme values of curva ture for cell sizes, and very little variation at higher cell sizes. The Kolmogorov-Smirnov (K-S) test was completed on each CFD for each terrain attribute. The rationale for this was to determine if there was a significant difference in the cumulative distribution of the terrain derivati ves when derived from DEMs of increasingly coarse resolutions. These results are presented in tabular format and their significance discussed below. 4.1.1 Elevation Figure 8 depicts the change in elevation representa tion along the two-mile transect. Represented elevation values are observe d to change and variability is reduced when represented in increasingly coarse DEMs. In s ome areas, represented elevations increase. At high resolutions, including those bet ween 20 and 160 feet, the elevations are represented in a very similar way, and there is fin e variability represented. The 320 foot DEM however begins to show a large loss of detail i n terrain variability, although abrupt and large elevation changes can be identified. Ele vation values are not realistically represented in the 640 foot DEM however; far too mu ch detail is lost, and even abrupt
56 and large changes in elevation are not easily ident ified. In this case, the cell size is simply too large to represent the fine variability. Over a length of about two miles, only 11 cells of 640 feet are used to represent the terr ain, whereas 528 cells of 20 feet are used over the same distance. Elevation Variation by Cell Size over Transect300 320 340 360 380 400 4200 2000400060008000 10000Distance (ft)Elevation (ft) Elevation 20 Elevation 40 Elevation 80 Elevation 160 Elevation 320 Elevation 640 Figure 8: Elevation Variability along a Transect a t Different Cell Sizes As indicated in Table 4, there is little change in the computed statistics for elevation. The mean value stays very constant with values of approximately 369 feet across all cell sizes. Additionally, there is litt le observed variation in the standard deviation across cell sizes as values range from 59 .18 feet for the 20 foot cell size, to a lower value of 57.9 feet at a cell size of 640 feet Larger effects of cell size are noted for the minimum and maximum elevation values at each ce ll size. The minimum elevation
57 value ranges from 203.87 feet for the 20 foot DEM t o 211.41 feet for the 640 foot cell size. Maximum elevation values exhibited a slightl y higher range of values, with the maximum value computed to be 519.1 feet computed fr om the 20 foot DEM, and 502.78 feet computed from the 640 foot DEM. This generall y indicates that maximum elevation becomes slightly underestimated when represented by a larger cell size DEM, and the minimum elevation is slightly overestimated when de rived from a larger cell size DEM. This is expected, as resampling uses the exact same elevation values for a much smaller number of cells. The relationship between represen ted elevation values and cell size is shown in Figure 9; indicating that elevation is not very resolution dependent, as expected. The trends in observed minimum and maximum elevatio n values are global. If just a single cell representing a maximum or minimum eleva tion value is not selected in the resampling, this could explain the change in maximu m and minimum values. Table 4: Summary statistics for elevation Cell Size (ft) Minimum Maximum Range Mean Standard Deviation 20 203.874 518.098 314.224 369.760 59.184 40 204.088 517.437 313.349 369.750 59.182 80 206.469 516.465 309.996 369.821 59.117 160 207.052 509.508 302.456 369.586 58.944 320 208.689 507.451 298.763 369.939 58.620 640 211.412 502.783 291.371 369.324 57.897
58 Minumum and Maximum Elevation Values for each Cell Size100 150 200 250 300 350 400 450 500 550 6000 4080 120160200240280320360400440480520560600640680Cell Size (ft)Elevation (ft) Maximum Elevation Minimum Elevation Figure 9: Minimum and Maximum Computed Elevation b y Cell Size Figure 10 shows the CFD for the elevation values in the watershed. As indicated by the graph, there is virtually no change in the cumulati ve distribution for the parameter. The same shape is shown in the distribution for each ce ll size. The cumulative frequency of elevation values do not change with a change in cel l size.
59 CFD Elevation0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 200250300350400450500 Elevation (ft)Cumulative Frequency Cell Size 20 Cell Size 40 Cell Size 80 Cell Size 160 Cell Size 320 Cell Size 640 Figure 10: CFD for Elevation Table 5: Results of the Kolmogorov-Smirnov Test on Elevation D Max Occurred At Total n Significance Elevation 20 and 40 0.0002 418.659 5812588 0.999 20 and 80 0.001 381.568 4940048 0.992 20 and 160 0.0027 394.664 4721865 0.629 20 and 320 0.0042 442.630 4667344 0.893 20 and 640 0.0116 399.436 4653712 0.528 As observed through the descriptive statistics, an d through an examination of the CFDs, there is very little difference in the mean e levation values, or the CFD for elevation. This is confirmed by the results of the K-S test for elevation. The D-Max values are extremely small between the CFD for each cell size, although they do increase progressively as the cell sizes increase. This ind icates an increasing difference in the distributions, though not enough to be significant, as indicated by the significance values. The CFD difference for elevation across cell sizes is statistically insignificant.
60 4.1.2 Slope The slopes represented by different cell sizes chan ge in response to the cell size from which they are derived. The steepest slope id entified in the profile is a slope of about 23 degrees, at a distance of about 5800 feet from the start of the transect, as shown in Figure 11. All of the peak slopes are identifie d from the 20 foot DEM. As cell size increases though, there is much less slope variabil ity present represented on the surface. Between cell sizes 20 and 80, the slopes represente d consistently decrease, through they still appear to represent the overall picture of sl ope on the surface, as compared to that derived from the 20 foot DEM. At cell sizes of 320 and 640 feet though, there is a distinct Â“flatteningÂ” of the terrain, and there is very little variability shown. This reduced representation of slope in the landscape has many h ydrologic implications, including changes in erosion potential, and potential soil sa turation.
61 Slope Variation by Cell Size over Transect0 5 10 15 20 250 2000400060008000 10000Distance (ft)Slope (degrees) Cell Size 20 Cell Size 40 Cell Size 80 Cell Size 160 Cell Size 320 Cell Size 640 Figure 11: Slope Variability along a Transect at D ifferent Cell Sizes The effect of grid cell size on slope is very subst antial as shown in Table 6 and Figure 12. There is a very obvious change in the m aximum, mean, and standard deviation values, indicating that slope decreases w hen computed from a DEM with a large cell size and these same DEMs are not able to identify areas of steep slope successfully. Maximum computed values of slope ran ge from a high of 56.91 degrees when computed from a 20 foot DEM to a low of 9.19 d egrees for the 640 DEM. Mean and standard deviation values also fall steadily fr om a high of 3.56 degrees to 2.98 degrees, and 3.03 degrees to 2.19 degrees respectiv ely. This relationship can also be readily understood from Figure 12. Mean slope fall s consistently from small cell sizes to
62 large ones, however, the maximum computed slope val ue falls sharply from small cell size to large ones. Table 6: Summary statistics for slope Cell Size (ft) Minimum Maximum Range Mean Standard Deviation 20 0.000 56.911 56.911 3.560 3.026 40 0.000 37.614 37.614 3.332 2.704 80 0.000 27.012 27.012 2.976 2.194 160 0.001 17.656 17.655 2.495 1.631 320 0.021 9.188 9.167 1.865 1.096 640 0.004 5.288 5.284 1.277 0.708 Slope Values for each Cell Size0 10 20 30 40 50 600 4080 120160200240280320360400440480520560600640680Cell Size (ft)Slope (degrees) Maximum Slope Mean Slope Figure 12: Maximum and Mean Computed Slope by Cell Size In order to clearly show the trend of mean slope va lues, the mean slope values for each cell size are shown separately in Figure 13, along with a log-trend line, formula and R2 value for the relationship.
63 Mean Slope Values for each Cell Sizey = -0.6717Ln(x) + 5.761 R2 = 0.97310 1 2 3 40 4080 120160200240280320360400440480520560600640680Cell Size (ft)Slope (degrees) Mean Slope Log. (Mean Slope) Figure 13: Mean Computed Slope by Cell Size with L og Trend Line Figure 13 shows a very strong log relationship betw een cell size and mean slope, as indicated through the very large R2 value of 0.9731. The CFD of slope is shown in Figure 14. Unlike el evation, the values for slope appear to be very sensitive to changes in the cell size from which they are derived. The arrow in Figure 14 indicates the leftward movement of the CFD to progressively lower values of slope as cell size increases. These curv es appear to be transformed versions of each other suggesting that the entire distribution, not just the mean, follows a relatively simple dependency on resolution. Slope is underest imated at large cell sizes, there is a smaller range of values of slope at large cell size s, and there are distinctly fewer high values of slope computed from larger cell sizes.
64 CFD Slope0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 05101520 Slope (degrees)Cumulative Frequency Cell Size 20 Cell Size 40 Cell Size 80 Cell Size 160 Cell Size 320 Cell Size 640 Figure 14: CFD for Slope Table 7: Results of the Kolmogorov-Smirnov Test on Slope Slope D Max Occurred At Total n Significance 20 and 40 0.0289 4.202062 5812588 0.001 20 and 80 0.0809 4.195947 4940048 0.001 20 and 160 0.177 3.600683 4721865 0.001 20 and 320 0.3405 2.769717 4667344 0.001 20 and 640 0.5294 2.112354 4653712 0.001 For the slope attribute, the D-Max statistic incre ases steadily as the cell size does, indicating an increasing difference in the CFD, and a greater maximum difference as cell size increases.
65 4.1.3 Overall Curvature Figure 15 shows the curvature values over the tran sect. As also evidenced through an examination of the descriptive statistic al values for the curvature parameters over the whole watershed, the curvature of the surf ace is greatly reduced, and the values of curvature become limited as cell size increases. When computed from a 20 foot DEM, there are many peaks and valleys in the curvature g raph, indicating relatively large values of negative and positive curvature, indicating a hi gh rate of change in slope for many locations over the transect. As cell sizes increas e, this variability quickly falls off, and values for the curvature parameters begin to stay v ery close to the value of zero; increasingly so as the cell size increases. At the cell size of 40 feet, there is limited change in the curvature and very little curvature e vident when computed from the 80 foot DEM. The graphed lines for the curvature parameter s appear to be indistinguishable from zero when computed from DEMs at 320 and 640 fe et.
66 Curvature Variation by Cell Size over Transect-2.5 -2.0 -1.5 -1.0 -0.5 0.0 0.5 1.0 1.5 2.0 2.50 2000400060008000 10000Distance (ft)Curvature (m m -2 ) Cell Size 20 Cell Size 40 Cell Size 40 Cell Size 160 Cell Size 320 Cell Size 640 Figure 15: Curvature Variability along a Transect at Different Cell Sizes Overall curvature values decrease greatly when com puted using large cell sizes. Descriptive statistics for curvature are given in T able 8. The mean curvature remains constant across all cell size ranges. This is expe cted, as globally the negative curvature values will balance with the positive values. The outer bound of minimum and maximum value for curvature as computed from each DEM becom e very limited in value. When computed from a 20 foot DEM, the minimum curvature is -21.82 m m-2 (upwardly concave) and the maximum curvature is 40.07 m m-2 (upwardly convex), indicating a terrain surface with extensive curvature. However, when computed from a 640 foot DEM, the minimum computed value is -0.035 m m-2 and the maximum value is 0.035 m
67 m-2 indicating a terrain surface with virtually no curv ature. Values for the standard deviation fall steadily as cell size increases, ind icating a reduction in the spread of different values, and the reduction in curvature va riability. Figure 16 shows the convergence of the curvature values as cell size in creases. While there is a large gap between minimum and maximum values at the 20 foot c ell size, this gap disappears as cell size increases, showing very little curvature representation when computed from 80 foot or larger DEMs. Table 8: Summary statistics for curvature Cell Size (ft) Minimum Maximum Range Mean Standard Deviation 20 -21.815 40.074 61.888 0.000001 0.331 40 -6.783 17.066 23.849 0.000069 0.197 80 -1.011 1.019 2.030 -0.000056 0.085 160 -0.330 0.479 0.808 -0.000152 0.044 320 -0.111 0.174 0.284 -0.000055 0.021 640 -0.035 0.035 0.070 -0.000051 0.008 Curvature Values for each Cell Size-25 -15 -5 5 15 25 35 450 4080 120160200240280320360400440480520560600640680Cell Size (ft)Curvature (m m-2) Maximum Curvature Minimum Curvature Figure 16: Maximum and Minimum Computed Curvature by Cell Size
68 Figure 17 show the cumulative frequencies of the v alues for curvature. The distribution shows a very similar pattern, just as was observed through an examination of the descriptive statistics of the values, and the v alues along a two-mile transect of the watershed. At small cell sizes, there are relative ly large values of curvature. These become increasingly smaller, and limited in their v ariability as cell sizes increase. The curvature parameter is strongly underestimated as c ell size increases, with values tending to zero as the cell size increases. With larger ce ll sizes in the DEM, it becomes very difficult to identify areas with large concavity va lues, or convexity values, both of which are important to hydrologic studies and hydrologic modeling as they are strong indicators of areas of flow convergence, and flow divergence a nd dispersion on the terrain surface. The arrows in Figure 17 indicate the direction of t he convergence of the values at increasingly coarse cell sizes. CFD Curvature0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 -1.0-0.50.00.51.0 Curvature Cumulative Frequency Cell Size 20 Cell Size 40 Cell Size 80 Cell Size 160 Cell Size 320 Cell Size 640 Figure 17: CFD for Curvature
69 Table 9: Results of the Kolmogorov-Smirnov Test on the Curvature Curvature D Max Occurred At Total n Significance 20 and 40 0.1093 0.12697600 5812588 0.001 20 and 80 0.2518 0.08297729 4940048 0.001 20 and 160 0.3354 0.06097412 4721865 0.001 20 and 320 0.4119 0.03697968 4667344 0.001 20 and 640 0.4736 0.01872253 4653712 0.001 The same trend is observed for the curvature attri bute as well, which is consistent with earlier observations the behavior of the attri bute through an examination of the elementary statistics and the graph of the CFD. As the cell size is increased, there is an increase in the maximum distance between the freque ncy distributions. 4.1.4 Plan Curvature Figure 18 shows the plan curvature values over the transect. As also evidenced through an examination of the descriptive statistic al values for the curvature parameters over the whole watershed, the plan curvature of the surface is greatly reduced, and the values of plan curvature become limited as cell siz e increases. When computed from a 20 foot DEM, there are many peaks and valleys in th e curvature graph, indicating relatively large values of negative and positive cu rvature, indicating a high rate of change in slope for many locations over the transect. As cell sizes increase, this variability quickly falls off, and values for the plan curvatur e parameters begin to stay very close to the value of zero; increasingly so as the cell size increases. At the cell size of 40 feet, there is limited change in the plan curvature and v ery little plan curvature evident when computed from the 80 foot DEM. The graphed lines f or the plan curvature parameter appear to be indistinguishable from zero when compu ted from DEMs at 320 and 640 feet.
70 Plan Curvature Variation by Cell Size over Transect-1.5 -1 -0.5 0 0.5 1 1.50 2000400060008000 10000Distance (ft)Plan Curvature (m m-2) Cell Size 20 Cell Size 40 Cell Size 80 Cell Size 160 Cell Size 320 Cell Size 640 Figure 18: Plan Curvature Variability along a Tran sect at Different Cell Sizes Plan curvature, a component of overall curvature s hows the same general pattern as seen for overall curvature. As shown in Table 1 0, as the cell size increases, there is an observed decrease in plan curvature represented in the terrain surface. As expected, the mean values for plan curvature remain steady with v alues very close to zero. The range of values as indicated by the standard deviation al so become very limited as cell size increases. Minimum and maximum values at the 20 fo ot cell size are smaller than seen for overall curvature, at -13.04 m m-2 (upwardly concave in the plan direction) and 25.82 m m-2 (upwardly convex in the plan direction) respectivel y, with values again near zero when computed from the 640 foot DEM. Effectively, the larger DEM cell sizes strongly limit the plan curvature on the surface. As with o verall curvature, the plan curvature
71 values converge when computed from larger cell size s, as shown in Figure 19. Above an 80 foot cell size, there is virtually no curvature represented in the DEM whatsoever. Table 10: Summary statistics for plan curvature Cell Size (ft) Minimum Maximum Range Mean Standard Deviation 20 -13.039 25.819 38.858 -0.000184 0.169 40 -3.932 8.401 12.333 0.000359 0.101 80 -0.715 0.761 1.477 0.000753 0.046 160 -0.310 0.297 0.607 0.000590 0.025 320 -0.077 0.091 0.167 0.000329 0.012 640 -0.018 0.022 0.040 0.000102 0.005 Plan Curvature Values for each Cell Size-15 -10 -5 0 5 10 15 20 25 300 4080 120160200240280320360400440480520560600640680Cell Size (ft)Plan Curvature (m m-2) Maximum Plan Curvature Minimum Plan Curvature Figure 19: Maximum and Minimum Computed Plan Curva ture by Cell Size Figure 20 show the cumulative frequencies of the v alues for plan curvature. The distribution shows a very similar pattern, just as was observed through an examination of the descriptive statistics of the values, and the v alues along a two-mile transect of the
72 watershed. At small cell sizes, there are relative ly large values of plan curvature. These become increasingly smaller, and limited in their v ariability as cell sizes increase. The curvature parameter is strongly underestimated as c ell size increases, with values tending to zero as the cell size increases. With larger ce ll sizes in the DEM, it becomes very difficult to identify areas with large concavity va lues, or convexity values, both of which are important to hydrologic studies and hydrologic modeling as they are strong indicators of areas of flow convergence, and flow divergence a nd dispersion on the terrain surface. The arrows in Figure 20 indicate the direction of t he convergence of the values at increasingly coarse cell sizes. CFD Plan Curvature0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 -1.0-0.50.00.51.0 Plan Curvature (m m-2)Cumulative Frequency Cell Size 20 Cell Size 40 Cell Size 80 Cell Size 160 Cell Size 320 Cell Size 640 Figure 20: CFD for Plan Curvature
73 Table 11: Results of the Kolmogorov-Smirnov Test o n the Plan Curvature Plan Curvature D Max Occurred At Total n Significan ce 20 and 40 0.1069 0.06827042 5812588 0.001 20 and 80 0.2348 0.04299930 4940048 0.001 20 and 160 0.3122 0.03305892 4721865 0.001 20 and 320 0.3906 0.02142961 4667344 0.001 20 and 640 0.4587 0.01077697 4653712 0.001 The same trend is observed for the plan curvature attribute as well, which is consistent with earlier observations the behavior o f the attribute through an examination of the elementary statistics and the graph of the C FD. As the cell size is increased, there is an increase in the maximum distance between the frequency distributions. 4.1.5 Profile Curvature Figure 21 shows the profile curvature values over the transect. As also evidenced through an examination of the descriptive statistic s values for the profile curvature parameters over the whole watershed, the profile cu rvature of the surface is greatly reduced, and the values of profile curvature become limited as cell size increases. When computed from a 20 foot DEM, there are many peaks a nd valleys in the profile curvature graph, indicating relatively large values of negati ve and positive curvature, indicating a high rate of change in slope for many locations ove r the transect. As cell sizes increase, this variability quickly falls off, and values for the profile curvature parameter begins to stay very close to the value of zero; increasingly so as the cell size increases. At the cell size of 40 feet, there is limited change in the pro file curvature and very little profile curvature evident when computed from the 80 foot DE M. The graphed lines for the
74 curvature parameters appear to be indistinguishable from zero when computed from DEMs at 320 and 640 feet. Profile Curvature Variation by Cell Size over Trans ect-2 -1.5 -1 -0.5 0 0.5 1 1.50 2000400060008000 10000Distance (ft)Profile Curvature (m m-2) Cell Size 20 Cell Size 40 Cell Size 80 Cell Size 160 Cell Size 320 Cell Size 640 Figure 21: Profile Curvature Variability along a T ransect at Different Cell Sizes Profile curvature, a component of overall curvatur e shows the same general pattern as seen for overall and plan curvature, as is expected. As the cell size increases, there is an observed decrease in computed profile c urvature represented in the terrain surface. As expected, the mean values for profile curvature remain steady with values very close to zero. The range of values as indicat ed by the standard deviation also become very limited as cell size increases. Minimu m and maximum values at the 20 foot cell size are smaller than seen for overall curvatu re, at -18.47 m m-2 (upwardly concave in
75 the profile direction) and 12.36 m m-2 (upwardly convex in the profile direction) respectively, with values again near zero when comp uted from the 640 foot DEM. Effectively, the larger DEM cell sizes strongly lim it the profile curvature on the surface. As with overall curvature, the profile curvature va lues converge when computed from larger cell sizes, as shown in Figure 22. Above an 80 foot cell size, there is virtually no curvature represented in the DEM whatsoever. Table 12: Summary statistics for profile curvature Cell Size (ft) Minimum Maximum Range Mean Standard Deviation 20 -18.471 12.363 30.834 -0.000185 0.216 40 -8.666 4.201 12.866 0.000290 0.125 80 -0.632 0.696 1.328 0.000809 0.054 160 -0.296 0.237 0.533 0.000742 0.027 320 -0.091 0.078 0.169 0.000385 0.012 640 -0.024 0.021 0.046 0.000153 0.005 Profile Curvature Values for each Cell Size-20 -15 -10 -5 0 5 10 150 4080 120160200240280320360400440480520560600640680Cell Size (ft)Profile Curvature (m m-2) Maximum Profile Curvature Minimum Profile Curvature Figure 22: Maximum and Minimum Computed Profile Cu rvature by Cell Size
76 Figure 23 show the cumulative frequencies of the v alues for profile curvature. The distribution shows a very similar pattern, just as was observed through an examination of the descriptive statistics of the va lues, and the values along a two-mile transect of the watershed. At small cell sizes, th ere are relatively large values of curvature. These become increasingly smaller, and limited in their variability as cell sizes increase. The profile curvature parameter is strongly underestimated as cell size increases, with values tending to zero as the cell size increases. With larger cell sizes in the DEM, it becomes very difficult to identify area s with large concavity values, or convexity values, both of which are important to hy drologic studies and hydrologic modeling as they are strong indicators of areas of flow convergence, and flow divergence and dispersion on the terrain surface. The arrows in Figure 23 indicate the direction of the convergence of the values at increasingly coars e cell sizes. CFD Profile Curvature0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 -1-0.500.51 Profile Curvature (m m-2)Cumulative Frequency Cell Size 20 Cell Size 40 Cell Size 80 Cell Size 160 Cell Size 320 Cell Size 640 Figure 23: CFD for Profile Curvature
77 Table 13: Results of the Kolmogorov-Smirnov Test o n the Profile Curvature Profile Curvature D Max Occurred At Total n Signifi cance 20 and 40 0.1071 -0.07232885 5812588 0.001 20 and 80 0.2483 -0.04712323 4940048 0.001 20 and 160 0.3359 -0.03424690 4721865 0.001 20 and 320 0.4131 -0.02057295 4667344 0.001 20 and 640 0.4724 -0.01124588 4653712 0.001 The same trend is observed for the profile curvatu re attribute as well, which is consistent with earlier observations the behavior o f the attribute through an examination of the elementary statistics and the graph of the C FD. As the cell size is increased, there is an increase in the maximum distance between the frequency distributions. 4.1.6 Stream Power Index The profile graph of the stream power index is sho wn in Figure 24. The results are most interesting when viewed in the context of the related terrain derivatives, especially elevation and slope. Beginning at about distance 2600 along the transect, there is a deep valley, as indicated in Figure 8, and is an area characterized by steep slopes as shown in Figure 11. In Figure 24, this same area i s shown to have a large value for the stream power index, at the 640 foot cell size. Thi s is also evident, at the smaller cell sizes, at locations of other smaller valleys at app roximately 6120 feet and 9000 feet along the transect. Generally, there seems to be more pe aks associated with the stream power index graphs at smaller cell sizes, most likely due to the fine variability of slopes and specific catchment area at these smaller cell sizes The spikes are generally not observed in the graph when stream power index is computed fr om larger cell sizes.
78 Stream Power Index Variation by cell Size over Tran sect0 50 100 150 200 250 300 350 4000 2000400060008000 10000Distance (ft)Stream Power Index Cell Size 20 Cell Size 40 Cell Size 80 Cell Size 160 Cell Size 320 Cell Size 640 Figure 24: Stream Power Index Variability along a Transect at Different Cell Sizes The resolution dependencies of the values for the stream power index are completely different from those noted for the other terrain derivatives, as reported in Table 14. The minimum values for the stream power index are constant at or near zero. The maximum values however, reflect the inverse of what is observed for both the mean and standard deviation values. For instance, the m aximum computed value for the stream power index is 4412.62 from the 20 foot DEM, with a steady fall in the computed values to a maximum value of 1300.89 observed from the 640 foot DEM. This behavior is the inverse of the observation for the mean stream powe r index values, which have a steady increase in value, from a low of 7.25 for the 20 fo ot DEM to a mean value of 29.81 for the 640 foot DEM. The standard deviation values in crease with an increase in cell size as
79 well, indicating a greater spread of value when com puted from a DEM with a larger cell size. The rise in mean values for the stream power index is not unexpected, as the values are proportional to the increased specific catchmen t area values computed from increasingly coarse DEMs. This steady increase in the mean value of stream power index with coarser DEMs is shown in Figure 25. While the difference in value is small between the 20 foot and 40 foot DEMs, and the 320 f oot and 640 foot DEMs, there is a much larger observed increase between the 80 foot, 160 foot, and 320 foot DEMs. The maximum stream power index value declines sharply f or the 320 and 640 foot DEMs. This correlates with the decrease in maximum slope observed in Table 6. The very low maximum slope values at lower resolutions could exp lain the large drop in maximum stream power index values as slope is directly rela ted to the stream power index. Table 14: Summary statistics for stream power inde x Cell Size (ft) Minimum Maximum Range Mean Standard Deviation 20 0.000 4412.630 4412.630 7.253640 31.771 40 0.000 2411.780 2411.780 9.601110 34.054 80 0.000 3821.980 3821.980 16.030400 59.481 160 0.001 3649.560 3649.560 23.195600 80.383 320 0.100 2060.890 2060.790 28.669100 78.388 640 0.045 1300.890 1300.850 29.805600 62.200
80 Mean Stream Power Index Values for each Cell Size5 10 15 20 25 300 4080 120160200240280320360400440480520560600640680Cell Size (ft)Stream Power Index Mean Stream Power Index Figure 25: Mean Computed Stream Power Index by Cel l Size As shown in Figure 25 above, there is a clear incr ease in the mean stream power index when computed from DEMs of increasingly coars e resolution. As noted in Table 14, there is a small increase in mean stream power index between the 20 foot and 40 foot cell size. Between the 40 foot and 80 foot cell si zes, and the 80 foot and 160 foot cell sizes, there is a much larger increase in mean stre am power index. As predicted, the stream power index increases as cell size increases This is a function of the increase in the specific catchment area in DEMs with large cell sizes. Above 320 feet, there is a very slight increase in the stream power index, as it se ems to reach a sill. Perhaps this is indicative that 409,600 square feet is near the lar gest specific catchment area, represented by a single 640 foot cell. Figure 26 shows the CFDs for the stream power inde x as computed from DEMs of each cell size. There is a general parallel shi ft to the right for the values, moving from
81 a small cell size to a large cell size, showing pro gressively higher values for the stream power index computed from DEMs with large cell size s. These curves confirm that the entire distribution is sensitive to cell size, not just the mean values. The frequency distribution becomes smoother as cell size increase s, not showing the large number of small values of stream power index as the frequency distribution for the values derived from the 20 foot DEM show. Generally, at each cell size, about 70% of the values appear to be under 20, with a much smaller percentage at a ll cell sizes having values above 20. The curves also shift at about this point for all c ell sizes, from a generally right-shift as cell size increases, to more of an upward shift; in dicating a smaller number of the very highest values. CFD Stream Power Index0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 0102030405060 Stream Power IndexCumulative Frequency Cell Size 20 Cell Size 40 Cell Size 80 Cell Size 160 Cell Size 320 Cell Size 640 Figure 26: CFD for Plan Stream Power Index
82 Table 15: Results of the Kolmogorov-Smirnov Test o n the Stream Power Index Stream Power Index D Max Occurred At Total n Significance 20 and 40 0.1453 0.7664566 5759502 0.001 20 and 80 0.3049 1.3066520 4897068 0.001 20 and 160 0.4615 1.9904480 4680470 0.001 20 and 320 0.5938 3.1079860 4626344 0.001 20 and 640 0.6743 4.0177690 4612825 0.001 An observation of the results of the D-Max statist ic for the stream power index indicate a trend similar to that seen for the other terrain attributes, where the D-Max value increases as cell size does. At each cell size, th e result of the difference is significant, meaning that the CFD of values is significantly dif ferent when computed from DEMs of increasingly coarse resolution. 4.1.7 Wetness Index Figure 27 shows the variation in the computed wetne ss index across the transect at various cell sizes. Generally, the surface is m odeled to be Â‘wetterÂ’ at larger cell sizes, and less prone to saturation at smaller cell sizes. The graph line from the 640 foot DEM is generally higher and much less variable than the graph line of the wetness index from smaller cell sizes. There is much more variability in the wetness index at the smaller cell sizes, but values are lower overall with most of th e values along the transect being smaller at the small cell sizes that the same locat ion represented by a larger cell size. This is consistent with the values observed in the descriptive statistics, where it was shown that the mean value of the wetness index incr eases at the watershed level as cell size does.
83 Wetness Index Variation by Cell Size over Transect4 6 8 10 12 14 16 180 2000400060008000 10000Distance (ft)Wetness Index Cell Size 20 Cell Size 40 Cell Size 80 Cell Size 160 Cell Size 320 Cell Size 640 Figure 27: Wetness Index Variability along a Trans ect at Different Cell Sizes Wetness index values vary quite strongly when compu ted from DEMs of increasingly coarse resolutions. There is a steady increase in both the mean and minimum values with larger cell sizes, as shown in Table 16. The minimum values range from 3.5 for the 20 foot DEM, to 9.96 for the 640 f oot DEM. There is also a steady increase in mean values from 8.7 for the 20 foot DE M to 12.41 for the 640 foot DEM. Although the mean values do not reflect this increa se directly, the maximum values behave in a different manner increasing steadily wh en computed from the 20 foot through the 160 foot DEM, and then falling slightly when co mputed from the 320 and 640 foot DEMs, to values seen at the finer resolutions. Gen erally though the wetness index values increase for the watershed with increasingly coarse resolutions, indicating the potential
84 for increased soil saturation, or the increased lik elihood of prediction of saturation, when using coarse DEMs to compute the wetness index. Th e relationship for the wetness index values are shown graphically in Figure 28. Interes tingly, the increase in minimum and mean wetness index values with an increase in cell size is nearly parallel, increasing most likely in a logarithmic fashion. Table 16: Summary statistics for wetness index Cell Size (ft) Minimum Maximum Range Mean Standard Deviation 20 3.494 19.301 15.807 8.700970 1.511 40 5.088 21.774 16.687 9.264280 1.485 80 6.242 23.706 17.464 9.949130 1.573 160 7.402 24.246 16.845 10.629800 1.625 320 8.693 19.755 11.062 11.465500 1.607 640 9.957 19.185 9.227 12.410000 1.450 Wetness Index Values for each Cell Size0 5 10 15 20 250 4080 120160200240280320360400440480520560600640680Cell Size (ft)Wetness Index Mean Wetness Index Min. Wetness Index Max. Wetness Index Figure 28: Maximum, Minimum, and Mean Computed Wet ness Index by Cell Size
85 In order to clearly show the trend of mean wetness index values, the mean wetness index values for each cell size are shown separately in F igure 29, along with a log-trend line, formula and R2 value for the relationship. Mean Wetness Index Values for each Cell Sizey = 1.0647Ln(x) + 5.3688 R2 = 0.9920 5 10 15 20 250 4080 120160200240280320360400440480520560600640680Cell Size (ft)Wetness Index Mean Wetness Index Log. (Mean Wetness Index) Figure 29: Mean Computed Wetness Index by Cell Siz e with Log Trend Line Figure 29 shows a very strong log relationship betw een cell size and mean wetness index, as indicated through the very large R2 value of 0.992. This logarithmic relationship is very similar to the form and strength as observed f or slope in Figure 13. The CFDs for the wetness index indicate a generall y parallel shift to the right as cell sizes increase, as shown in Figure 30. These higher values indicate an increased potential or prediction of soil saturation. The po tential for soil saturation increases as the cell size from which it is computed does. There ar e very few areas of very high wetness index, and no abrupt shifts in the distribution. G enerally though, these distributions show
86 that the surface is predicted to be wetter when rep resented by DEMs with increasingly coarse resolutions. CFD Wetness Index0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.03456789 101112131415161718Wetness IndexCumulative Frequency Cell Size 20 Cell Size 40 Cell Size 80 Cell Size 160 Cell Size 320 Cell Size 640 Figure 30: CFD for Wetness Index Table 17: Results of the Kolmogorov-Smirnov Test o n the Wetness Index Wetness Index D Max Occurred At Total n Significanc e 20 and 40 0.1937 8.309067 5810087 0.001 20 and 80 0.4019 8.678735 4937846 0.001 20 and 160 0.5818 9.068766 4719673 0.001 20 and 320 0.7528 9.712436 4665152 0.001 20 and 640 0.8723 10.573270 4651520 0.001 The results for the wetness index parameter are co mparable to that of the stream power index. While the shapes of the distribution are different, there is similarity in the way that the D-Max statistic increases as cell size does. When compared to the CFD for the 20 foot cell size, the distribution for each su ccessive distribution is significantly different, and the D-Max statistic increases in val ue. This indicates an increase in the
87 difference between the distributions as cell size i ncreases. For all terrain attributes except elevation, the results indicate that there is a sta tistically significant difference between the CFDs of values when computed from large DEM cell si zes, as compared to those computed from the high resolution 20 foot DEM.
88 4.2 Stream Analysis 4.2.1 Stream Network Comparison Figures 31 through 36 show a visual overall compari son of the stream networks as mapped through an orthophoto (Figure 31) and as der ived from a DEM for successively coarse DEM resolutions (Figures 32 through 36). Mapped StreamsOrder 1 2 3 4 5 6 01234 0.5 Miles Figure 31: Orthophoto Derived Stream Network
89 DEM Derived 20 ftOrder 1 2 3 4 5 6 01234 0.5 Miles Figure 32: Streams Derived from 20 foot DEM
90 DEM Derived 40 ftOrder 1 2 3 4 5 6 01234 0.5 Miles Figure 33: Streams Derived from 40 foot DEM
91 DEM Derived 80 ftOrder 1 2 3 4 5 6 01234 0.5 Miles Figure 34: Streams Derived from 80 foot DEM
92 DEM Derived 160 ftOrder 1 2 3 4 5 6 01234 0.5 Miles Figure 35: Streams Derived from 160 foot DEM
93 DEM Derived 320 ftOrder 1 2 3 4 5 6 01234 0.5 Miles Figure 36: Streams Derived from 320 foot DEM Some difference in the stream networks is apparent when viewed at a small scale, however, viewing the network at a large scale is ne cessary to observe the effect of DEM resolution upon the derivation of a stream networks from DEMs of increasingly coarse resolutions. Figures 37 to 42 below show, for a sa mple of the study area, how the stream network differs when derived from a DEM of increasi ngly coarse resolution, with the same contributing area threshold, as compared to a photogrammetrically mapped stream network. Figure 37 shows the streams derived from a 20 foot DEM compared to photogrammetrically mapped streams. Visually, the stream networks match very closely, with similar meanders, both following the valleys e vident on the surface very well. There
94 is some difference in stream length. In some cases the photogrammetrically mapped streams are slightly longer than the derived stream s, and in other cases the derived streams are longer. When determining the contribut ing area threshold to use to initiate a stream network, one was chosen that resulted in a s tream network that most closely resembled the mapped stream network in terms of the channel initiation points. This threshold was found by experimenting with a range o f threshold values. The derived stream network, even from a high resolution DEM wil l not be consistent in all streams when a single stream initiation threshold is used. DEM Derived 20' Orthophoto Streams Figure 37: Streams Derived from a 20 foot Cell Siz e and Photogrammetrically Mapped Streams
95 As the streams are derived from DEMs of coarser ce ll sizes, detail is lost in the stream network, as shown in Figure 38. This Figure shows the stream network derived from a 40 foot DEM, and the mapped stream network. Although similar in morphology to the mapped stream network, and the stream networ k derived from the 20 foot DEM, it is obvious that some detail is lost in the stream n etwork. DEM Derived 40' Orthophoto Streams Figure 38: Streams Derived from a 40 foot Cell Siz e and Photogrammetrically Mapped Streams Figure 39 shows the relationship between the mappe d stream network, and that derived from a DEM with a cell size of 80 feet. At this resolution, much detail in the stream network is lost. The streams generally foll ow the depressions in the DEM, and follow the general morphology of the mapped stream network. Although stream lengths
96 appear to be similar between the derived streams an d the mapped streams, the fine detail in the network observed at higher resolutions is go ne. There are many straight stream segments, which are very different than the mapped stream network. The network becomes more linear and straight at this resolution resulting in only a vague representation of the natural stream network. DEM Derived 80' Orthophoto Streams Figure 39: Streams Derived from a 80 foot Cell Siz e and Photogrammetrically Mapped Streams The generalization of the stream network increases when derived from a 160 foot DEM, as compared to both the photogrammetrically ma pped stream network, and those derived from higher resolution DEMs. As shown in F igure 40, the stream network only vaguely resembles the mapped stream network. The d erived stream network consists of
97 straight segments with very little detail. Any fin e meander or detail that is present in the mapped network is not present in the derived networ k, a function of the fact that the meander cannot be represented with such a large cel l size, and of the generalization that would naturally occur in a stream network derived f rom such a coarse DEM. At this cell size, the derived network crosses through peaks in the surface where a stream would not be expected to flow, and is not representative of t he natural stream network. DEM Derived 160' Orthophoto Streams Figure 40: Streams Derived from a 160 foot Cell Si ze and Photogrammetrically Mapped Streams The final cell size for which a stream network was derived is 320 feet. Streams were not derived from the 640 foot DEM, because at this cell size the initiation threshold area is less than the area of one cell at 640 feet. The result of the derivation is shown in
98 Figure 41. At this cell size, any detail that was present in the mapped stream network, or that derived from the high resolution DEM is lost. Streams are very general, and the lengths are very underestimated, as the streams are straight lines in the derived network, where they meander in nature. The derived stream n etwork follows very generally the morphology of the mapped network, but because of th e large cell size, the stream network is very far from its natural location. Thi s is a function of the cell size used to derive the network; the natural network simply cann ot be represented adequately at this cell size. Many small streams are not represented in this derived network either. Quite possibly, the natural streams that are not represen ted are shorter than the 320 foot cell size that they are being derived from.
99 DEM Derived 320' Orthophoto Streams Figure 41: Streams Derived from a 320 foot Cell Si ze and Photogrammetrically Mapped Streams The final Figure in this series, Figure 42 shows a ll of the derived stream networks together with the photogrammetrically mapped stream s. It is obvious that there is much difference between each of the derived networks and the photogrammetrically derived networks. The differences in the networks increase as the cell size from which they are derived does. At all cell sizes though, there are some small stream segments that are not represented in the derived network whatsoever. Whe n derived from large cell sizes, the location of the stream networks tend to not represe nt the true location of the network, nor do they represent the true morphology of the natura l stream network. Streams can only be located in very general terms when derived from DEMs with large cell sizes.
100 DEM Derived 20' DEM Derived 40' DEM Derived 160' DEM Derived 80' DEM Derived 320' Orthophoto Streams Figure 42: All Streams Derived from a DEM and Phot ogrammetrically Mapped Streams 4.2.2 Stream Network Statistical Analysis To evaluate the accuracy of the streams derived dir ectly from a 20 foot DEM compared to streams mapped through photogrammetric methods, several techniques were used. The Kappa Index of Agreement (KIA) was used to investigate the stream agreement, by order, between these two stream datas ets. To complete this task, the photogrammetrically mapped streams were ordered usi ng the RivEx stream ordering tool, and then rasterized to a 20 foot cell size while be ing snapped to the 20 foot DEM so that the cells would be aligned in the same grid system as the DEM derived stream network.
101 The KIA is able to test for raster independence, an d quantify the level of raster agreement. Table 18 gives the cell count for the interception of all streams from the 20 foot DEM, irrespective of order. The overwhelming numbe r of cells in each raster did not represent stream networks at all. Direct compariso ns at this level are not useful, there are simply too many cells that are not streams overall. However, an examination of the KIA is useful for comparison, as the KIA corrects for c hance, comparing the cell agreement against the agreement that might be expected as cha nce. For this reason, the KIA can be thought of as the chance-corrected proportion of ag reement (Chuang, 2001). Table 18: Cell Agreement Count Matrix for Boolean Stream Network Derived Streams Mapped 0 1 Total Streams 0 17715936 129085 17845021 1 93325 61654 154979 Total 17809261 190739 18000000 Table 19 gives the raw cell agreement counts betwee n the rasterized photogrammetrically mapped streams, and the streams derived from the 20 foot DEM, by order. The numbers from Tables 18 and 19 are then used to derive the K IA. The KIA for Tables 18 and 19 are both different, but both are meaningful. The v alues in Table 18 are a measure of agreement without considering order. The values in Table 19 are a measure of agreement considering order. The agreement values in Table 1 9 are used to determine the KIA for each order of stream. These values are also used t o calculate a third value of the KIA, one that does not explicitly determine the KIA by o rder, but for the stream network as a whole considering order.
102 Table 19: Cell Agreement Count Matrix for 20' DEM Derived vs. Mapped Stream Network Mapped Stream Network 0 1 2 3 4 5 6 Total Derived 0 17715936 50696 20492 11117 4258 2910 151 17805560 Stream 1 71046 19416 2983 425 205 228 231 94534 Network 2 31371 7920 9692 1282 132 151 68 50616 3 14288 1207 4163 4822 729 81 172 25462 4 5683 143 143 1717 2016 12 129 9843 5 3476 231 34 26 418 1456 0 5641 6 3221 80 65 24 1 2 1297 4690 Total 17845021 79693 37572 19413 7759 4840 2048 17996346 Table 20 gives the KIA for the photogrammetrically mapped streams and the streams derived from the 20 foot DEM. Results are given for the overall stream network irrespective of stream order, the overall stream ne twork considering order, and for each stream order individually. Higher numbers for the KIA indicate a greater strength of agreement. Generally, the KIA results indicate a l ow level of agreement between the two stream networks. The agreement for the overall net work at 0.35 is higher than any other value besides the sixth-order stream agreement. Th is seems likely, as the KIA in this case was based strictly on whether a cell was a str eam in both rasters, not on whether the corresponding order was captured correctly. The lo w value of 0.238 for the first order network is expected as well. As observed previousl y, because the DEM derived streams are based on a contributing area threshold, the loc ation of initiation of a first-order stream is often incorrect. Sometimes the point of initiat ion of a first-order stream in a DEM derived network will correspond directly with that of a Â“blue lineÂ” stream in a mapped network. Other times, the location for initiation of a first-order stream will be different, as shown in Figure 37. This low value for the KIA reflects this fact.
103 DEM Derived 20 ftOrder 1 2 3 4 5 6 Mapped StreamsOrder 1 2 3 4 5 6 02805608401,120 140 Feet Figure 43: Result of Error in Initiating a First-O rder Stream, Cascade Effect Figure 43 shows the effect on stream order when a f irst-order stream is initiated in the wrong location. The two short stream segments in t he inlaid box are classified as firstorder streams when derived from the 20 foot DEM. H owever, these stream segments are not present in the mapped stream network. Instead, the second-order derived stream segment begins at approximately the same location a s the first-order mapped stream. This effect will cascade down the stream network, m isclassifying the downstream segments. As a comparison, the first-order streams in the lower right corner of Figure 43 show the initiation of the first-order derived stre am as being very close to the initiation point of the mapped first-order stream.
104 It is logical that the KIA generally increases as the stream order value increases. While the location of initiation of a first-order s tream is often incorrect, there is greater chance that the location of higher-order streams wi ll agree. Visually, there is less variability in the stream networks as the stream or der increases, and this is reflected in the KIA values in Table 20. The KIA for the overall ne twork considering order is 0.289. This value is lower than the KIA value of 0.35 that was observed for the stream network not considering order, as expected. This again ind icates that there is some misclassification of stream order. When the stream networks are considered as Boolean networks, there is more agreement, as there is no s tream order classification to consider. Table 20: Kappa Index of Agreement for 20 foot DEM Kappa Index of Agreement (KIA) Statistic Stream Location Agreement for 20' DEM Overall KIA without considering order 0.350 Overall KIA with consideration of order 0.289 Individual KIA First Order 0.238 Individual KIA Second Order 0.232 Individual KIA Third Order 0.274 Individual KIA Fourth Order 0.288 Individual KIA Fifth Order 0.313 Individual KIA Sixth Order 0.434 4.2.3 Comparison of Vector Stream Locations To further investigate the agreement between the ma pped stream network, and those derived from DEMs of varying resolutions, an analysis of the vector stream networks was completed. This involved buffering th e mapped stream network by the
105 width of one pixel (20 feet) and intersecting the f ive different stream networks with this buffer. The result is an interception table that g ives the percentage of the network intercepted by the buffer for each cell size, and f or each order. The results are presented in Table 21. Table 21: Results of Intersection of Mapped Stream Buffer and Derived Stream Networks DEM Cell Size (ft) 20 40 80 160 320 Overall Interception 67.31% 65.54% 64.24% 60.06% 59.76% First Order 32.75% 32.27% 32.15% 30.28% 30.11% Second Order 28.94% 28.61% 29.54% 27.93% 27.64% Third Order 22.75% 24.64% 23.12% 22.91% 26.71% Fourth Order 22.59% 27.01% 20.73% 22.20% 20.94% Fifth Order 33.03% 35.09% 29.96% 28.18% 29.02% Sixth Order 42.83% 41.21% 37.00% 36.04% 29.04% Overall, the length of the network intercepted decl ines steadily between the 20 foot and 320 foot cell sizes. The overall stream network is successfully intercepted 67.31% of its length, while this decreases to only 59.76% for the stream network generated with the 320 foot DEM. As is evident from these observation s, a stream network more closely resembling the mapped network is derived from the h igher resolution DEM. The results of the network interception by order a re not always consistent across all resolutions. For first-order streams, about 32 % of the stream network is intercepted for all resolutions, falling to about 30% for the 1 60 and 320 foot steam networks. The results for second-order streams are similar. Abou t 29% of the second-order stream network is intercepted for the 20, 40, and 80 foot DEMs, falling slightly by about 2% for the 160 and 320 foot cell sizes. For the third-ord er streams, although the percentage of the stream network intercepted is similar, between about 22% and 24% intercepted for the 20 to 160 foot cell size, the amount of interce ption is largest for the 320 foot cell size
106 at about 26%. This contrasts the first and secondorder stream networks that have their lowest interception values at the higher cell sizes The results for the fourth-order streams are inconsistent. For the 20, 40, 80, 160, and 320 foot cell sizes respectively, the results are 22.59%, 27.01%, 20.73%, 22.2%, and 20.9 4%. These results show no clear trend, but the amount of interception is generally low. The results for fifth-order streams do not show a clear trend either. The greatest amo unt of interception is observed at with the 20 foot and 40 foot cell size stream networks, at about 33% and 35% respectively. This represents a large improvement in interception over the fourth-order networks, and a much higher level of interception than observed for the 80, 160, and 320 foot stream networks. For the sixth-order streams, the interce ption rate is generally much higher, a result consistent with the results of the KIA discu ssed earlier. Sixth-order streams generated from the 20 foot DEM were intercepted cor rectly almost 43% of the time. The results for the 40 foot stream network were similar with about a 41% interception rate. The rate of interception decreased steadily for the larger cell sizes though, to a low of only 29% for the 320 foot stream network. Generall y, it appears that the higher resolution stream network more closely resembles th e photogrammetrically mapped stream network. It is important to note however, t hat the overall interception is much better than the interception observed by order. A likely reason for this is the use of a constant stream initiation threshold used to initia te streams. First-order streams will often start in the incorrect place, resulting in the firs t-order derived streams linking up with second-order photogrammetrically mapped streams or vice-versa. This error will cascade down the rest of the stream orders as seen in Figur e 43. As was observed with the KIA,
107 stream accuracy results are in part affected by str eam length issues, and the location of initiation of first-order streams, second-order str eams, etc. 4.2.4 Stream Morphology Statistics In addition to reporting the level of stream agreem ent between the derived stream network and the photogrammetrically mapped stream n etwork, the stream morphology characteristics provide a number of good metrics fo r comparison. The stream morphology statistics provide information about the stream network that is not apparent from an examination of the location of the stream n etworks alone. For instance, Table 22 gives the total count of the number of stream segme nts, by order, for each of the stream networks including the photogrammetrically mapped s tream network from Wake County. Results for the first-order stream networks reveal that there are many more first order streams in each of the derived networks other than the 320 foot network, than there is in the mapped stream network. For the 20 foot c ell size, there are 968 more first order streams present than in the mapped stream network, a significant increase. In general terms, there is a steady decrease in the number of first order streams derived as the cell size increases. This increase could be partly due again to the choice of a constant contributing area threshold used to initiate a stre am network from the DEM. Removal of short first-order streams from the derived stream n etwork could result in the number of first-order streams being closer to those in the ma pped network. As discussed previously, when a constant channel initiation threshold is use d, the initiation point of a channel can often be incorrect, and streams that could be simpl y ephemeral streams or gullies could be incorrectly identified as first-order streams. Â“PruningÂ” these short first-order streams
108 from the network could reduce the possibility of th is error. This consideration is more likely when the results for the second-order stream network counts are observed. Further, the channel initiation technique could be reconside red, for instance a variable threshold could be used instead of a constant stream initiati on threshold. As expected, there are far fewer second-order stre ams at each cell size than there are first-order streams. However, the results are much more tightly grouped for the second order streams. For instance, there are 1280 second order stream segments in the 20 foot stream network, whereas there are 900 obser ved in the photogrammetrically mapped stream network for a difference of only 380. The number of second-order stream segments decrease steadily for each increased cell size, with the number of second-order streams being underestimated at the 320 foot cell s ize by 238 stream segments. The number of stream segments are underestimated by the 320 foot stream network for every order, except for the fourth-order, where the resul ts are extremely close; within one stream segment of each other. Overall, the numbers of stream segments in the der ived networks tend to overestimate the actual number of segments at every order, with the exception of the 320 foot network, which underestimates. Also in genera l terms, as the cell size from which a stream network is derived increases, the trend is f or the number of stream segments by order to decrease. By increasing the cell size tha t is used to derive a stream network, the detail in that network begins to decrease, perhaps making the stream network less reliable. Importantly, higher DEM resolution does not provide better estimates of stream segment counts; however the results are very depend ent upon the channel initiation technique used.
109 Table 22: Count of Stream Segments by Order Cell Size (ft) First Order Second Order Third Order Fourth Order Fifth Order Sixth Order Total 20 2784 1280 779 329 218 156 5546 40 2722 1233 759 309 207 161 5391 80 2572 1172 628 306 215 159 5052 160 2645 1113 696 246 211 140 5051 320 1796 662 385 200 92 111 3246 Mapped Streams 1816 900 439 199 149 164 3667 The bifurcation ratio of a stream network gives th e relationship between the number of streams in one order to the number of str eams in the next successive order. In nature, bifurcation ratios of between 3 and 5 are c ommon (Summerfield, 1991). Table 23 summarizes the bifurcation ratios computed for each derived stream network and the photogrammetrically mapped stream network. The res ults for the bifurcation ratios by order don't reveal consistent increases or decrease s, but indicate a less reliable estimate of stream network properties. None of the derived str eam networks represent the bifurcation ratios in the photogrammetrically mappe d streams very well. Importantly, higher DEM resolution does not provide better estim ates of the bifurcation ratios, much the same as was observed with raw stream segment co unts. The results are very dependent however upon the channel initiation techn ique used.
110 Table 23: Bifurcation Ratios Cell Size (ft) RB1 RB2 RB3 RB4 RB5 20 2.18 1.64 2.37 1.51 1.40 40 2.21 1.62 2.46 1.49 1.29 80 2.19 1.87 2.05 1.42 1.35 160 2.38 1.60 2.83 1.17 1.51 320 2.71 1.72 1.93 2.17 0.83 Mapped Streams 2.02 2.05 2.21 1.34 0.91 Table 24 presents the actual length of streams in each order. As observed in the table, as DEM cell size increases, the overall stre am length decreases. This trend applies to all stream orders except for first-order streams Length of first-order streams remain roughly the same as cell size increases. These str eam lengths are used to compute the stream length ratio which is another metric used to make comparisons of streams in one order to streams in another order. In this case, t he total length of streams in one order was compared to the total length of streams in the next successive order. The results are presented in Table 25. The results for stream leng th ratios show fairly substantial changes with resolution, with a substantial decreas e overall. Table 24: Length of Streams by Order (in feet) Cell Size (ft) First Order Second Order Third Order Fourth Order Fifth Order Sixth Order Total 20 1572989 851686 427970 165051 95911 77950 3191557 40 1559209 806728 417413 150066 85718 78599 3097733 80 1580573 777858 346017 153285 92066 78343 3028143 160 1731869 709981 359460 123860 83998 70610 3079777 320 1618614 551407 263184 124641 54033 69062 2680942 Mapped Streams 1266056 606765 313032 125860 77837 93927 2483476
111 Table 25: Stream Length Ratios Cell Size (ft) RL1 RL2 RL3 RL4 RL5 20 0.54 0.50 0.39 0.58 0.81 40 0.52 0.52 0.36 0.57 0.92 80 0.49 0.44 0.44 0.60 0.85 160 0.41 0.51 0.34 0.68 0.84 320 0.34 0.48 0.47 0.43 1.28 Mapped Streams 0.48 0.52 0.40 0.62 1.21 The results in Table 25 indicate that that higher r esolution DEMs do not produce better estimates of stream length order, but the results a re very sensitive to the channel initiation technique used to derive them.
112 Chapter 5 Conclusions Beven and Moore (1992) noted that digital terrain m odeling is not new in the field of hydrology, and in fact the characterization of t he topography must be important in hydrology. However, what is different now is that the era of GIS, simulation visualization software, and workstation computing p ower has arrived and is rapidly maturing. New, higher resolution elevation data fr om new collection methods and sensors such as LIDAR has provided the opportunity to explore hydrologic and other physical processes at scales not previously explore d. It also raises a question as to the advantage (if any) of high resolution elevation dat a in the study of hydrologic and terrain processes. It is often preferable to allow necessi ty to drive technology, rather than have technology influence our needs. Are high resolutio n elevation data sets a case where technology has led the change? This thesis research had one goal; to determine th e effect of DEM resolution upon hydrologically significant terrain attributes, incl uding stream networks. It proposed to determine this effect through two objectives. The first objective was to determine the effect of DEM resolution on hydrologically signific ant terrain derivatives including slope, plan curvature, profile curvature, overall c urvature, topographic wetness index, and the stream power index. The second objective was t o determine how the resolution of a DEM affects the morphology and attributes of the st reams derived from it. Measures to
113 investigate this include the Kappa Index of Agreeme nt (KIA), the length of derived streams intercepted by a buffered photogrammetrical ly mapped stream network, and the comparison of several morphological characteristics of the derived stream network, to that of the photogrammetrically mapped stream netwo rk. To facilitate this research, a high resolution DEM data set derived from LIDAR data was acquired for a small North Carolina waters hed. This DEM was resampled to successively coarser resolutions resulting in DEMs with 20, 40, 80, 160, 320, and 640 foot cell sizes. Terrain attributes with significa nce to hydrology and hydrologic modeling were then computed and compared between the resolut ions. DEM cell size was observed to have a significant e ffect upon the value of the computed attributes. As hypothesized, as the cell size from which slope was computed increased, the computed value of slope decreased. The amount of decrease was large; a maximum slope of 56.9 degrees was observed in the w atershed when represented by a 20 foot DEM, while the maximum slope fell to only 5.3 degrees when represented by a 640 foot DEM. The mean slope values fell progressively as well, following a logarithmic curve. A logarithmic trend line added to the chart of mean slope graphed with DEM cell size resulted in a R2 value of 0.9731. Slope has implications for overl and and subsurface flow velocity and runoff rate and affects precipita tion infiltration, among other surface processes. Values for overall, plan, and profile surface curv ature were hypothesized to become limited in their spread, and the mean values of the curvature attributes were expected to decrease with an increase in cell size from which they are computed. This phenomenon was observed though examination of the s ummary statistics and cumulative
114 frequency diagrams for the attribute. Generally, w hen computed from a DEM of 80 foot cell size or greater, there is very little curvatur e represented in the surface and the results have little meaning for terrain analysis or hydrolo gic modeling. The curvature attribute has implications for flow acceleration, the erosion and deposition rate, converging and diverging flow, soil water content, and soil charac teristics. The stream power index and the topographic wetness index were topographic attributes included in the research for a number of reasons. These included the fact that they are both compound terrain attributes, encompas sing a first derivative attribute, and because of their importance to hydrology and hydrol ogic modeling. The stream power index was hypothesized to increase as cell size inc reased, as this attribute will increase in proportion to an increase in the specific catchment area, which increases with an increase in DEM cell size used to derive it. The stream pow er index increased steadily with an increase in cell size, from a value of 7.25 from a 20 foot cell size to a high of 29.8 when computed from a 604 foot cell size DEM. This attri bute has importance to hydrology, as it is able to predict areas of net erosion and net deposition when combined with a curvature parameter indicating areas of flow accele ration or deceleration. As previous studies (Zhang and Montgomery, 1994) h ave concluded and as was hypothesized in this research, the topographic wetn ess index increases in value, and causes a surface to be modeled as Â‘wetterÂ’ when com puted from DEMs of increasing resolution. The topographic wetness index describe s the spatial distribution and extent of zones of for runoff generation as a function of ups lope contributing area, soil transmissivity, and slope. It is also a key compon ent of the TOPMODEL (Beven and Kirkby, 1979) hydrologic modeling framework.
115 The second objective of this research was to deter mine the effect of DEM resolution on the stream networks derived from DEMs To test this, streams were delineated using the same contributing area thresho ld from a DEM with cell sizes of 20, 40, 80, 160, and 320 feet. These were then compare d to a photogrammetrically mapped stream network. Visual and quantitative comparison s of resulting stream locations were completed, and certain morphological attributes wer e compared to determine stream agreement. The Kappa Index of Agreement (KIA) was also used t o determine the amount of agreement for the 20 foot stream network. Agreemen t values were generally low with an overall agreement of 0.35, and values by order rang ing from 0.23 to 0.43. The accuracy of each derived stream network was determined by co mparing it to the photogrammetrically mapped stream network. The acc uracy of the stream network was expected to decrease when derived from DEMs with la rge cell sizes. This was confirmed through an observation of the percentage of the str eam network intercepted by the buffered photogrammetrically mapped stream network. Interception rates generally fell overall and by stream order when derived from DEMs of increasingly coarse resolutions. Overall interception was 67.3% for the 20 foot stre am network, falling to only a 59.8% interception rate for the 320 foot stream network. In addition to stream location agreement, stream m orphology metrics were investigated. Bifurcation ratios were calculated f or each of the derived stream networks and the photogrammetrically mapped stream network. The bifurcation ratios were not expected to decrease between one stream order and t he next successive one as the cell size increased, but remain fairly constant. The re sults for the bifurcation ratios by order
116 don't reveal consistent increases or decreases, but indicate a less reliable estimate of stream network properties. None of the derived str eam networks represent the bifurcation ratios in the photogrammetrically mappe d streams very well. For all stream orders other than first-order strea ms, as hypothesized, the stream length decreased as cell size used to derive the st ream increased. This indicates less meander in the stream and a straightening of the st ream network. This was also confirmed through a visual inspection of the stream network derived at each resolution. First-order stream length can be somewhat misleadin g, as the length of a first-order stream will change with the channel initiation tech nique used to initiate it, making these lengths often incorrect. This research indicates that stream networks deriv ed from larger cell size DEMs are less reliable than when derived from high-resol ution DEMs. Through the observation of stream morphology and accuracy characteristics, it is noted that steam networks derived from high-resolution DEMs are able to more closely reflect a natural stream network. Derived stream networks more closely refl ect the mapped natural stream network as stream-order increases. This was observ ed through the values of the stream interception analysis, with derived sixth-order str eams matching the natural network the closest. The same results were obtained through th e KIA. The derived sixth-order streams showed the most agreement with the mapped n atural stream network. The physical location of derived streams becomes l ess accurate when derived from DEMs with increasingly coarse resolutions, but there is little difference in the overall stream network statistics. For example, me andering is represented poorly using coarse DEMs, but this is true across all orders. F or this reason, metrics like the
117 bifurcation ratio are not necessarily less accurate when computed from low resolution DEMs as compared to high resolution DEMs. These re sults indicate that when available, a high-resolution DEM would be preferable for strea m network modeling. Clearly, it is difficult to reliably derive the lo cation of initiation of a first-order stream. This is further complicated by the use of a constant contributing area threshold across a whole watershed to determine the beginning of the stream network. This was a major reason why the stream networks derived from t he LIDAR DEMs were not closer to the mapped streams. The correct length and count o f first-order streams is difficult to determine, and most likely often incorrect. There are many factors that can influence the valu es of terrain attributes computed from DEMs. Of these, cell size is only one factor. LIDAR data collection specifications, interpolation techniques, filtering techniques, flo w routing algorithm, contributing area threshold, and error are a number of factors. All of these factors need to be considered when interpreting topographic attributes derived fr om DEMs and subsequently using the results to model hydrological processes. Sensor technology continues to improve, and comput ing power continues to allow the mass consumption of very large, high resolution data sets. Based on the results of this research, using the highest resolution elevation da ta available is beneficial to terrain analysis and stream network analysis. Technology m ight be pushing the possibilities forward, but it is up to the terrain analyst or hyd rologist to quantify and discover the benefit of these new data sources for their area of interest. Perhaps the ultimate goal of this research is to establish relationships that qu antify the effects of DEM resolution upon hydrologically relevant terrain derivatives, which can then be considered when
118 processing DEMs from various resolutions for the pu rpose of parameterizing hydrologic models. The future of this research could include using a DEM with even higher resolution. Now emerging are even higher resolutio n LIDAR data, including data with a cell size of less than 10 feet. Terrestrial based LIDAR sensors are bound to decrease the cost of LIDAR elevation acquisition and increase th e availability and use of LIDAR derived elevation data sets. Using alternatives to the contributing area thresh old to initiate streams is a research opportunity provided by LIDAR DEMs. As ob served, the channel initiation technique used to initiate a stream network can hav e a profound impact on the results of subsequent analysis, especially in terms of stream network morphology. By developing more adaptive techniques that are more flexible in terms on stream initiation, perhaps more accurate stream networks can be derived from D EMs where there is no mapped stream network, and more reliable stream morphology statistics can be generated. The development of empirical relationships to desc ribe the scale-dependency of terrain attributes is a future direction of this re search. The results of this study indicate that empirical relationships can be derived between DEM cell size and the values of terrain attributes derived from DEMs. What is nece ssary to carry this work further is to research different locations at different scales to confirm what these relationships look like. This research could be conducted across a wid e range of terrain types, i.e. coastal, plains, low, moderate, and high relief areas and lo cations representing areas of net erosion like the Appalachian Mountains, and areas o f net deposition, like coastal Florida. Instead of a single watershed, a number of watershe ds with different relief and area could
119 be compared. While this research was conducted wit hin a sixth-order watershed, would the same results be observed in a watershed of grea ter order? Further, resolution independent terrain derivative s need to be developed. There are few, if any, of these developed at the present time. Much research effort has been directed at quantifying the effect of DEM cell size upon terrain derivatives. The body of knowledge surrounding this is increasing rapidly, e specially in this era of higher resolution data and increased computing power. Dev eloping resolution independent metrics and derivatives is a logical next step in t his process.
120 References Cited Adler, D. 2001. IBM DB2 Spatial Extender Â– Spatial data within the RDBMS. In Proceedings of the 27th VLDB Conference, Roma, Italy Aronoff, S. 1989. Geographic Information Systems: A Management Persp ective WDL Publications. Atkinson, P. 2002. Transactions in GIS 6(1): 1-4. Band, L. E. 1986. Topographic partition of waters heds with digital elevation models. Water Resources Research 22(1): 15-24. Beven, K. J and Kirkby, M. J. 1979. A physically based, variable contributing area model of basin hydrology. Hydrological Sciences Bulletin 24: 43-69. Beven, K. J., and Moore, I. D., editors. 1992. Terrain Analysis and Distributed Modelling in Hydrology John Wiley & Sons. Brasington, J., and Richards, K. 1998. Interactio ns between model predictions, parameters, and DTM scales for TOPMODEL. Computers & Geosciences 24(4): 299-314. Burrough, P.A. 1986. Principles of Geographical Information Systems for Land Resources Assessment Oxford University Press, New York. Burrough P. A. and McDonnell, R. A. 1998. Principles of Geographical Information Systems New York, Oxford University Press. Carter, J.R. 1988. Digital representations of top ographic surfaces. Photogrammetric Engineering and Remote Sensing 54: 1577-80. Chairat, S., Delleur, J.W., 1993. Integrating a physically based hydrological model w ith GRASS In Kovar, K., Nachtnebel, H. P., editors. HydroGI S 93: Application of Geographical Information Systems in Hydrology and Water Resources International Association of Hydrological Sciences, Vienna.
121 Chuang, Jen-Hsiang. 2001. Agreement between Catego rical Measurements: Kappa Statistics. Retrieved 24 Mar 2006, from Kappa Stat istics. Website: http://www.dmi.columbia.edu/homepages/chuangj/kappa / Costa-Cabral, M.C. and Burges, S. J. 1994. Digital elevation model networks (DEMON): A model of flow over hillslopes for computation of contributing and dispersal areas. Water Resources Research 30: 1681-92. Creed, I. E. and Band, L. E. 1998. Exploring simila rity in the export behavior of nitrate-H from forested catchments: A mechanistic modeling a pproach, Water Resources Research 34(11): 3079-3094. Desmet P. J. and Govers, G. 1996. Comparison of rou ting algorithms for Digital Elevation Models and their implications for predict ing ephemeral gullies. International Journal of Geographical Information S ystems 10: 311-31. Doan, J. 2000. Hydrologic Model of the Buffalo Ba you Using GIS. In Maidment D R and Djokic D (Editors) Hydrologic and Hydraulic Modeling Support with Geographic Information Systems Redlands, CA, ESRI Press: 113-143. Duke, G., Kienzle, S. W., Johnson, D. L. and Byrne, J. M. 2003. Improving overland flow routing by incorporating ancillary road data into Digital Elevation Models. Journal of Spatial Hydrology 3(2). ESRI. 2005. ArcGIS Documentation. Fairfield, J. and Leymarie, P. 1991. Drainage netwo rks from grid Digital Elevation Models. Water Resources Research 27: 709-17. Garbrecht, J. and Martz, L. W. 2000. Digital elevat ion model issues in water resources modeling. In Maidment, D. R. and Djokic, D., edit ors. Hydrologic and Hydraulic Modeling Support with Geographic Information Syste ms ESRI Press: 1-28. Goodchild, M. 2004. GIScience, Geography, Form, a nd Process. Annals of the Association of American Geographers 94(4): 709-714. Hellweger, F. L. 1997. AGREE Â– DEM Surface Recond itioning System. Retrieved 24 Mar 2006, from AGREE Â– DEM Surface Reconditioning System. Website: http://www.ce.utexas.edu/prof/maidment/gishydro/fe rdi/research/agree/ agree.html Horn, B.K.P. 1981. Hill-shading and the reflectance map. Proceedings of the IEEE 19 (1), 14-47.
122 Hutchinson, M. F. 1989. A new procedure for gridd ing elevation and stream line data with automatic removal of spurious pits. Journal of Hydrology 106: 211-232. Hutchinson, M. F. and Gallant, J. C. 2000. Digita l Elevation Models and Representation of Terrain Shape. In Terrain Analysis: Principles and Applications John Wiley & Sons, Inc. Jenson, S. 1992. Applications of Hydrologic Infor mation Automatically Extracted from Digital Elevation Models. In Beven, K. J., and Mo ore, I. D. editors. Terrain Analysis and Distributed Modelling in Hydrology John Wiley & Sons. Jensen, S. and Domingue, J. O. 1988. Extracting to pographic structure from digital elevation model data for geographic information sy stem analysis. Photogrammetric Engineering and Remote Sensing 54: 1593-1600. Jones, K. H. 1998. A comparison of two approaches to ranking algorithms used to compute hill slopes. GeoInformatica 2:3, 235-256. Kavaya, M. 1999. LIDAR Tutorial. Retrieved 24 Ma r 2006 from LIDAR Tutorial. Website: http://www.ghcc.msfc.nasa.gov/sparcle/spa rcle_tutorial.html/. Kenny, F. and Matthews, B. 2005. A methodology fo r aligning raster flow direction data with photogrammetrically mapped hydrology. Computers & Geosciences 31: 768-779. Kienzle, S. 2004. The Effect of DEM Raster Resolu tion on First Order, Second Order, and Compound Terrain Derivatives. Transactions in GIS 8(1): 83-111. Lindsay J. B. 2005. The Terrain Analysis System: A tool for hydro-geomorphic applications. Hydrological Processes 19(5): 1123-1130. Lindsay, J. B. and Creed, I. F. 2005. Removal of artifact depressions from digital elevation models: towards a minimum impact approac h. Hydrological Processes 19(16): 3113-3126. Longley, P.A., Goodchild, M. F., Maguire, D. J., an d Rhind, D. W. 2001. Geographic Information Systems and Science. New York: Wiley. Marceau, D. J. and Hay, G. J. 1999. Remote sensing contributions to the scale issue. Canadian Journal of Remote Sensing 25(4): 357-366. Melville, J. K. and Martz, L. W. 2004. A Comparis on of Data Sources for Manual and Automated Hydrographical Network Delineation. Can adian Water Resources Journal 29(4): 267-282.
123 Moore, I. D., Grayson, R. B. and Ladson, A. R. 1991 Digital Terrain Modeling: a review of hydrological, geomorphological, and biol ogical applications. Hydrological Processes 5: 3-30. National Research Council. 1991, Opportunities in the Hydrological Sciences National Academy Press, Washington, D.C. NOAA. LIDAR Data Handler: A NOAA Coastal Services Center ArcView or ArcMap Extension. Retrieved 24 March 2006 from Topograph ic Change Mapping Â– LIDAR Data Handler. Website: htttp://www.csc.noaa.gov/crs/tcm/lidar_handler.htm l North Carolina Floodplain Mapping Program. 2003. L IDAR and Digital Elevation Fact Sheet, January 2003. O'Callaghan, J. F. and Mark, D. M. 1984. The extrac tion of drainage networks from digital elevation data. Computer Vision, Graphics and Image Processing 28: 32844. Peuker, T. K. and Douglas, D. H. 1975. Detection of surface-specific points by local parallel processing of discrete terrain elevation data. Computer Graphics and Image Processing 4: 375-387. Plaster, R. L. Geospatial Solutions. February 200 2. Quinn P. F., Beven, K., Chevallier, P. and Planchon O. 1991. The prediction of hillslope flow paths for distributed hydrological modeling u sing digital terrain models. H ydrological Processes 5: 59-79. Saulnier, G., Obled, C., and Beven, K. 1997. Analyt ical compensation between DTM grid resolution and effective values of saturated hydraulic conductivity within the TOPMODEL framework. Hydrological Processes 11: 1331 46 Saunders, W. K. and Maidment, D. R. 1996. A GIS Assessment of Nonpoint Source Pollution in the San Antonio-Nueces Coastal Basin Center for Research in Water Resources Online Report 96-1, University of Texas, Austin, TX. Schneider, B. 2001. Phenomenon-based Specificatio n of the Digital Representation of Terrain Surfaces. Transactions in GIS 5(1): 39-52 Strahler, A. N. 1957. Quantitative Analysis of Wate rshed Geomorphology. Transactions of the American Geophysical Union 8(6): 913-920. Sole, A., and Valanzano, A. In V. P. Singh and M. Fiorentine(editors) 1996. Geographical Information Systems in Hydrology. 175-194.
124 Summerfield, M. 1991. Global Geomorphology London: Longman Scientific and Technical. Tarboton, D. G. 1991. On the Extraction of Channel Networks from Digital Elevation Data. Hydrological Processes 5: 81-100. Tarboton, D. G. 1996. Fractal river networks, Hort onÂ’s laws and Tokunaga cyclicity. Journal of Hydrology 187: 105-117. Tarboton, D. G. 1997. A New Method for the Determin ation of Flow Directions and Contributing Areas in Grid Digital Elevation Model s. Water Resources Research 33(2): 309-319. Tarboton, D. G. 2004. Terrain Analysis Using Digi tal Elevation Models. Retrieved 24 Mar 2006 From TauDEM Terrain Analysis Using Digita l Elevation Models. Website: http://hydrology.neng.usu.edu/taudem/ Wilson, J. P. and Gallant, J. C. 2000. Digital Te rrain Analysis. In Terrain Analysis: Principles and Applications John Wiley & Sons, Inc. Wolock, D. and McCabe, G. 2000. Differences in to pographic characteristics computed from 100and 1000-m resolution digital elevation model data. Hydrological Processes 14: 987-1002. Wolock, D. and Price, C. 1994. Effects of digital elevation model map scale and data resolution on a topography-based watershed model. Water Resources Research 30(11): 3041-3052 Zhang, W. and Montgomery, D. R. 1994. Digital eleva tion model grid size, landscape representation, and hydrologic simulations. Water Resources Research 30: 1019 28. Zhou and Liu. 2002. Error assessment of grid-based flow routing algorithms used in hydrological models. International Journal of Geographical Information S cience 16(8): 819-842.
126 Appendix A: Geoprocessing Scripts Python script used to create the flow direction and flow accumulation rasters # -------------------------------------------------------------------------# Create_fdr_and_fac_Rasters.py # -------------------------------------------------------------------------# Import system modules import sys, string, os, win32com.client # Create the Geoprocessor object gp = win32com.client.Dispatch("esriGeoprocessing.Gp Dispatch.1") # Check out any necessary licenses gp.CheckOutExtension("spatial") # Load required toolboxes... gp.AddToolbox("C:/Program Files/ArcGIS/ArcToolbox/T oolboxes/Spatial Analyst Tools.tbx") # Local variables... fdr640 = "D:\\ThesisDataFolder_Streams\\fdr640" Output_drop_raster = "" fdr320 = "D:\\ThesisDataFolder_Streams\\fdr320" Output_drop_raster__2_ = "" fdr160 = "D:\\ThesisDataFolder_Streams\\fdr160" Output_drop_raster__3_ = "" fdr80 = "D:\\ThesisDataFolder_Streams\\fdr80" Output_drop_raster__4_ = "" fdr40 = "D:\\ThesisDataFolder_Streams\\fdr40" Output_drop_raster__5_ = "" filldem40 = "D:\\ThesisDataFolder_Streams\\filldem4 0" filldem80 = "D:\\ThesisDataFolder_Streams\\filldem8 0" filldem160 = "D:\\ThesisDataFolder_Streams\\filldem 160" filldem320 = "D:\\ThesisDataFolder_Streams\\filldem 320" filldem640 = "D:\\ThesisDataFolder_Streams\\filldem 640" fac40 = "D:\\ThesisDataFolder_Streams\\fac40" fac80 = "D:\\ThesisDataFolder_Streams\\fac80" fac160 = "D:\\ThesisDataFolder_Streams\\fac160"
127 Appendix A (continued) fac320 = "D:\\ThesisDataFolder_Streams\\fac320" fac640 = "D:\\ThesisDataFolder_Streams\\fac640" # Process: Flow Direction... gp.FlowDirection_sa(filldem640, fdr640, "NORMAL", O utput_drop_raster) # Process: Flow Direction (2)... gp.FlowDirection_sa(filldem320, fdr320, "NORMAL", O utput_drop_raster__2_) # Process: Flow Direction (3)... gp.FlowDirection_sa(filldem160, fdr160, "NORMAL", O utput_drop_raster__3_) # Process: Flow Direction (4)... gp.FlowDirection_sa(filldem80, fdr80, "NORMAL", Out put_drop_raster__4_) # Process: Flow Direction (5)... gp.FlowDirection_sa(filldem40, fdr40, "NORMAL", Out put_drop_raster__5_) # Process: Flow Accumulation... gp.FlowAccumulation_sa(fdr40, fac40, "") # Process: Flow Accumulation (2)... gp.FlowAccumulation_sa(fdr80, fac80, "") # Process: Flow Accumulation (3)... gp.FlowAccumulation_sa(fdr160, fac160, "") # Process: Flow Accumulation (4)... gp.FlowAccumulation_sa(fdr320, fac320, "") # Process: Flow Accumulation (5)... gp.FlowAccumulation_sa(fdr640, fac640, "")
128 Appendix A (continued) Python script used to import topographic wetness in dex and stream power index ASCII files from the Terrain Analysis System to Arc GIS # -------------------------------------------------------------------------# Import_WI_SPI_from_TAS.py # -------------------------------------------------------------------------# Import system modules import sys, string, os, win32com.client # Create the Geoprocessor object gp = win32com.client.Dispatch("esriGeoprocessing.Gp Dispatch.1") # Load required toolboxes... gp.AddToolbox("C:/Program Files/ArcGIS/ArcToolbox/T oolboxes/Conversion Tools.tbx") # Local variables... WI640 = "D:\\ThesisDataFolder_Derivatives\\wi640" RSP640 = "D:\\ThesisDataFolder_Derivatives\\rsp640" WI320 = "D:\\ThesisDataFolder_Derivatives\\wi320" RSP320 = "D:\\ThesisDataFolder_Derivatives\\rsp320" WI160 = "D:\\ThesisDataFolder_Derivatives\\wi160" RSP160 = "D:\\ThesisDataFolder_Derivatives\\rsp160" WI80 = "D:\\ThesisDataFolder_Derivatives\\wi80" RSP80 = "D:\\ThesisDataFolder_Derivatives\\rsp80" WI40 = "D:\\ThesisDataFolder_Derivatives\\wi40" RSP40 = "D:\\ThesisDataFolder_Derivatives\\rsp40" WI20 = "D:\\ThesisDataFolder_Derivatives\\wi20" RSP20 = "D:\\ThesisDataFolder_Derivatives\\rsp20" lidar20asc_RSP_asc = "D:\\ThesisDataFolder_Derivati ves\\TAS\\lidar20asc_RSP.asc" lidar20asc_WI_asc = "D:\\ThesisDataFolder_Derivativ es\\TAS\\lidar20asc_WI.asc" lidar40asc_RSP_asc = "D:\\ThesisDataFolder_Derivati ves\\TAS\\lidar40asc_RSP.asc" lidar40asc_WI_asc = "D:\\ThesisDataFolder_Derivativ es\\TAS\\lidar40asc_WI.asc" lidar80asc_RSP_asc = "D:\\ThesisDataFolder_Derivati ves\\TAS\\lidar80asc_RSP.asc" lidar80asc_WI_asc = "D:\\ThesisDataFolder_Derivativ es\\TAS\\lidar80asc_WI.asc" lidar160asc_RSP_asc = "D:\\ThesisDataFolder_Derivat ives\\TAS\\lidar160asc_RSP.asc" lidar160asc_WI_asc = "D:\\ThesisDataFolder_Derivati ves\\TAS\\lidar160asc_WI.asc" lidar320asc_RSP_asc = "D:\\ThesisDataFolder_Derivat ives\\TAS\\lidar320asc_RSP.asc" lidar320asc_WI_asc = "D:\\ThesisDataFolder_Derivati ves\\TAS\\lidar320asc_WI.asc" lidar640asc_RSP_asc = "D:\\ThesisDataFolder_Derivat ives\\TAS\\lidar640asc_RSP.asc" lidar640asc_WI_asc = "D:\\ThesisDataFolder_Derivati ves\\TAS\\lidar640asc_WI.asc"
129 Appendix A (continued) # Process: ASCII to Raster... gp.ASCIIToRaster_conversion(lidar640asc_WI_asc, WI6 40, "FLOAT") # Process: ASCII to Raster (2)... gp.ASCIIToRaster_conversion(lidar640asc_RSP_asc, RS P640, "FLOAT") # Process: ASCII to Raster (3)... gp.ASCIIToRaster_conversion(lidar320asc_WI_asc, WI3 20, "FLOAT") # Process: ASCII to Raster (4)... gp.ASCIIToRaster_conversion(lidar320asc_RSP_asc, RS P320, "FLOAT") # Process: ASCII to Raster (5)... gp.ASCIIToRaster_conversion(lidar160asc_WI_asc, WI1 60, "FLOAT") # Process: ASCII to Raster (6)... gp.ASCIIToRaster_conversion(lidar160asc_RSP_asc, RS P160, "FLOAT") # Process: ASCII to Raster (7)... gp.ASCIIToRaster_conversion(lidar80asc_WI_asc, WI80 "FLOAT") # Process: ASCII to Raster (8)... gp.ASCIIToRaster_conversion(lidar80asc_RSP_asc, RSP 80, "FLOAT") # Process: ASCII to Raster (9)... gp.ASCIIToRaster_conversion(lidar40asc_WI_asc, WI40 "FLOAT") # Process: ASCII to Raster (10)... gp.ASCIIToRaster_conversion(lidar40asc_RSP_asc, RSP 40, "FLOAT") # Process: ASCII to Raster (11)... gp.ASCIIToRaster_conversion(lidar20asc_WI_asc, WI20 "FLOAT") # Process: ASCII to Raster (12)... gp.ASCIIToRaster_conversion(lidar20asc_RSP_asc, RSP 20, "FLOAT")
130 Appendix A (continued) Python script used to resample the original LIDAR D EM and generate the topographic attributes of slope, and the three curv ature attributes, and convert the results to ASCII format for subsequent import to th e Terrain Analysis System # -------------------------------------------------------------------------# Model_Resample_Gen_Derivatives.py # -------------------------------------------------------------------------# Import system modules import sys, string, os, win32com.client # Create the Geoprocessor object gp = win32com.client.Dispatch("esriGeoprocessing.Gp Dispatch.1") # Check out any necessary licenses gp.CheckOutExtension("spatial") # Load required toolboxes... gp.AddToolbox("C:/Program Files/ArcGIS/ArcToolbox/T oolboxes/Spatial Analyst Tools.tbx") gp.AddToolbox("C:/Program Files/ArcGIS/ArcToolbox/T oolboxes/Conversion Tools.tbx") gp.AddToolbox("C:/Program Files/ArcGIS/ArcToolbox/T oolboxes/Data Management Tools.tbx") # Local variables... lidar640ras = "D:\\ThesisDataFolder_Derivatives\\li dar640ras" lidar_20_fl = "lidar_20_fl" lidar40ras = "D:\\ThesisDataFolder_Derivatives\\lid ar40ras" lidar320ras = "D:\\ThesisDataFolder_Derivatives\\li dar320ras" lidar80ras = "D:\\ThesisDataFolder_Derivatives\\lid ar80ras" lidar160ras = "D:\\ThesisDataFolder_Derivatives\\li dar160ras" slope40 = "D:\\ThesisDataFolder_Derivatives\\slope4 0" slope80 = "D:\\ThesisDataFolder_Derivatives\\slope8 0" slope160 = "D:\\ThesisDataFolder_Derivatives\\slope 160" slope320 = "D:\\ThesisDataFolder_Derivatives\\slope 320" slope640 = "D:\\ThesisDataFolder_Derivatives\\slope 640" Slope20 = "D:\\ThesisDataFolder_Derivatives\\slope2 0" curvature640 = "D:\\ThesisDataFolder_Derivatives\\c urvature640" profile640 = "D:\\ThesisDataFolder_Derivatives\\pro file640" plan640 = "D:\\ThesisDataFolder_Derivatives\\plan64 0" curvature320 = "D:\\ThesisDataFolder_Derivatives\\c urvature320" profile320 = "D:\\ThesisDataFolder_Derivatives\\pro file320"
131 Appendix A (continued) plan320 = "D:\\ThesisDataFolder_Derivatives\\plan32 0" curvature160 = "D:\\ThesisDataFolder_Derivatives\\c urvature160" profile160 = "D:\\ThesisDataFolder_Derivatives\\pro file160" plan160 = "D:\\ThesisDataFolder_Derivatives\\plan16 0" curvature20 = "D:\\ThesisDataFolder_Derivatives\\cu rvature20" profile20 = "D:\\ThesisDataFolder_Derivatives\\prof ile20" plan20 = "D:\\ThesisDataFolder_Derivatives\\plan20" curvature40 = "D:\\ThesisDataFolder_Derivatives\\cu rvature40" profile40 = "D:\\ThesisDataFolder_Derivatives\\prof ile40" plan40 = "D:\\ThesisDataFolder_Derivatives\\plan40" curvature80 = "D:\\ThesisDataFolder_Derivatives\\cu rvature80" profile80 = "D:\\ThesisDataFolder_Derivatives\\prof ile80" plan80 = "D:\\ThesisDataFolder_Derivatives\\plan80" slope40asc_asc = "D:\\ThesisDataFolder_Derivatives\ \slope40asc.asc" profile40asc_asc = "D:\\ThesisDataFolder_Derivative s\\profile40asc.asc" plan40asc_asc = "D:\\ThesisDataFolder_Derivatives\\ plan40asc.asc" curvature40asc_asc = "D:\\ThesisDataFolder_Derivati ves\\curvature40asc.asc" slope80asc_asc = "D:\\ThesisDataFolder_Derivatives\ \slope80asc.asc" profile80asc_asc = "D:\\ThesisDataFolder_Derivative s\\profile80asc.asc" plan80asc_asc = "D:\\ThesisDataFolder_Derivatives\\ plan80asc.asc" curvature80asc_asc = "D:\\ThesisDataFolder_Derivati ves\\curvature80asc.asc" slope20asc_asc = "D:\\ThesisDataFolder_Derivatives\ \slope20asc.asc" profile20asc_asc = "D:\\ThesisDataFolder_Derivative s\\profile20asc.asc" plan20asc_asc = "D:\\ThesisDataFolder_Derivatives\\ plan20asc.asc" curvature20asc_asc = "D:\\ThesisDataFolder_Derivati ves\\curvature20asc.asc" slope160asc_asc = "D:\\ThesisDataFolder_Derivatives \\slope160asc.asc" profile160asc_asc = "D:\\ThesisDataFolder_Derivativ es\\profile160asc.asc" plan160asc_asc = "D:\\ThesisDataFolder_Derivatives\ \plan160asc.asc" curvature160asc_asc = "D:\\ThesisDataFolder_Derivat ives\\curvature160asc.asc" slope320asc_asc = "D:\\ThesisDataFolder_Derivatives \\slope320asc.asc" profile320asc_asc = "D:\\ThesisDataFolder_Derivativ es\\profile320asc.asc" plan320asc_asc = "D:\\ThesisDataFolder_Derivatives\ \plan320asc.asc" curvature320asc_asc = "D:\\ThesisDataFolder_Derivat ives\\curvature320asc.asc" slope640asc_asc = "D:\\ThesisDataFolder_Derivatives \\slope640asc.asc" profile640asc_asc = "D:\\ThesisDataFolder_Derivativ es\\profile640asc.asc" plan640asc_asc = "D:\\ThesisDataFolder_Derivatives\ \plan640asc.asc" curvature640asc_asc = "D:\\ThesisDataFolder_Derivat ives\\curvature640asc.asc" lidar20asc_asc = "D:\\ThesisDataFolder_Derivatives\ \lidar20asc.asc" lidar80asc_asc = "D:\\ThesisDataFolder_Derivatives\ \lidar80asc.asc" lidar40asc_asc = "D:\\ThesisDataFolder_Derivatives\ \lidar40asc.asc" lidar640asc_asc = "D:\\ThesisDataFolder_Derivatives \\lidar640asc.asc" lidar320asc_asc = "D:\\ThesisDataFolder_Derivatives \\lidar320asc.asc" lidar160asc_asc = "D:\\ThesisDataFolder_Derivatives \\lidar160asc.asc"
132 Appendix A (continued) # Process: Resample (2)... tempEnvironment0 = gp.workspace gp.workspace = "D:\\ThesisDataFolder_Derivatives" gp.Resample_management(lidar_20_fl, lidar40ras, "40 ", "NEAREST") gp.workspace = tempEnvironment0 # Process: Slope 40... gp.Slope_sa(lidar40ras, slope40, "DEGREE", "1") # Process: Raster to ASCII... gp.RasterToASCII_conversion(slope40, slope40asc_asc ) # Process: Curvature_40... gp.Curvature_sa(lidar40ras, curvature40, "1", profi le40, plan40) # Process: Raster to ASCII (2)... gp.RasterToASCII_conversion(profile40, profile40asc _asc) # Process: Raster to ASCII (3)... gp.RasterToASCII_conversion(plan40, plan40asc_asc) # Process: Raster to ASCII (4)... gp.RasterToASCII_conversion(curvature40, curvature4 0asc_asc) # Process: Resample (4)... gp.Resample_management(lidar_20_fl, lidar80ras, "80 ", "NEAREST") # Process: Slope 80... gp.Slope_sa(lidar80ras, slope80, "DEGREE", "1") # Process: Raster to ASCII (5)... gp.RasterToASCII_conversion(slope80, slope80asc_asc ) # Process: Curvature_80... gp.Curvature_sa(lidar80ras, curvature80, "1", profi le80, plan80) # Process: Raster to ASCII (6)... gp.RasterToASCII_conversion(profile80, profile80asc _asc) # Process: Raster to ASCII (7)... gp.RasterToASCII_conversion(plan80, plan80asc_asc)
133 Appendix A (continued) # Process: Raster to ASCII (8)... gp.RasterToASCII_conversion(curvature80, curvature8 0asc_asc) # Process: Slope 20... gp.Slope_sa(lidar_20_fl, Slope20, "DEGREE", "1") # Process: Raster to ASCII (9)... gp.RasterToASCII_conversion(Slope20, slope20asc_asc ) # Process: Curvature_20... gp.Curvature_sa(lidar_20_fl, curvature20, "1", prof ile20, plan20) # Process: Raster to ASCII (10)... gp.RasterToASCII_conversion(profile20, profile20asc _asc) # Process: Raster to ASCII (11)... gp.RasterToASCII_conversion(plan20, plan20asc_asc) # Process: Raster to ASCII (12)... gp.RasterToASCII_conversion(curvature20, curvature2 0asc_asc) # Process: Resample (5)... gp.Resample_management(lidar_20_fl, lidar160ras, "1 60", "NEAREST") # Process: Slope 180... gp.Slope_sa(lidar160ras, slope160, "DEGREE", "1") # Process: Raster to ASCII (13)... gp.RasterToASCII_conversion(slope160, slope160asc_a sc) # Process: Curvature_160... gp.Curvature_sa(lidar160ras, curvature160, "1", pro file160, plan160) # Process: Raster to ASCII (14)... gp.RasterToASCII_conversion(profile160, profile160a sc_asc) # Process: Raster to ASCII (15)... gp.RasterToASCII_conversion(plan160, plan160asc_asc ) # Process: Raster to ASCII (16)... gp.RasterToASCII_conversion(curvature160, curvature 160asc_asc) # Process: Resample (3)...
134 Appendix A (continued) gp.Resample_management(lidar_20_fl, lidar320ras, "3 20", "NEAREST") # Process: Slope 320... gp.Slope_sa(lidar320ras, slope320, "DEGREE", "1") # Process: Raster to ASCII (17)... gp.RasterToASCII_conversion(slope320, slope320asc_a sc) # Process: Curvature_320... gp.Curvature_sa(lidar320ras, curvature320, "1", pro file320, plan320) # Process: Raster to ASCII (18)... gp.RasterToASCII_conversion(profile320, profile320a sc_asc) # Process: Raster to ASCII (19)... gp.RasterToASCII_conversion(plan320, plan320asc_asc ) # Process: Raster to ASCII (20)... gp.RasterToASCII_conversion(curvature320, curvature 320asc_asc) # Process: Resample... gp.Resample_management(lidar_20_fl, lidar640ras, "6 40", "NEAREST") # Process: Slope 640... gp.Slope_sa(lidar640ras, slope640, "DEGREE", "1") # Process: Raster to ASCII (21)... gp.RasterToASCII_conversion(slope640, slope640asc_a sc) # Process: Curvature_640... gp.Curvature_sa(lidar640ras, curvature640, "1", pro file640, plan640) # Process: Raster to ASCII (22)... gp.RasterToASCII_conversion(profile640, profile640a sc_asc) # Process: Raster to ASCII (23)... gp.RasterToASCII_conversion(plan640, plan640asc_asc ) # Process: Raster to ASCII (24)... gp.RasterToASCII_conversion(curvature640, curvature 640asc_asc) # Process: Raster to ASCII (25)... gp.RasterToASCII_conversion(lidar_20_fl, lidar20asc _asc)
135 Appendix A (continued) # Process: Raster to ASCII (26)... gp.RasterToASCII_conversion(lidar80ras, lidar80asc_ asc) # Process: Raster to ASCII (27)... gp.RasterToASCII_conversion(lidar40ras, lidar40asc_ asc) # Process: Raster to ASCII (28)... gp.RasterToASCII_conversion(lidar640ras, lidar640as c_asc) # Process: Raster to ASCII (29)... gp.RasterToASCII_conversion(lidar320ras, lidar320as c_asc) # Process: Raster to ASCII (30)... gp.RasterToASCII_conversion(lidar160ras, lidar160as c_asc)