Handling Gridded Data Topography and Projections GIS in

Handling Gridded Data: Topography and Projections GIS in Water Resources Fall 2014 by Ayse Kilic with materials from David G. Tarboton, Utah State University and from ESRI software

Learning Objectives • Calculation of slope on a raster using – Arc. GIS method based in finite differences – D 8 steepest single flow direction – D steepest outward slope on grid centered triangular facets • Map Projections – State Plane Coordinate System – UTM (Universal Transverse Mercator Coordinate System)

Spatial Surfaces used in Hydrology Elevation Surface — the ground surface elevation at each point -- Expressed as a Digital Elevation Model for Gridded Data

Types of Elevation Data available Data Spatial reference • GCS_WGS_1984 • Decimal degrees (Global Topography) • WGS 1984 GTOPO SRTM (Shuttle Radar Topography Mission) NED 30 (National Elevation Data) NED 10 Z Bit Depth units 16 -bit 30 arcsec signed/un m (1 km) signed Integer Pixel size • GCS_WGS_1984 • Decimal degrees • WGS 1984 90 m m 16 -bit signed Integer • GDC_North_America_1983 • Decimal degrees • NAD 1983 1 arcsec (30 m) m Float • GDC_North_America_1983 • Decimal degrees • NAD 1983 1/3 arcsec m (10 m) Float 3 ft Float • NAD 83_HARN_State. Plane_Oregon_North Lidar • Foot (DEM/DSM) • NAD 1983 HARN ft

Slope Handout http: //snr. unl. edu/kilic/giswr/2014/Slope. pdf Determine the length, slope and azimuth of the line AB.

3 -D detail of the Tongue river at the WY/Montana border from LIDAR from aircraft or from the ground can provide amazing detail on elevation, including individual tree heights and hydraulic channels Roberto Gutierrez University of Texas at Austin

Topographic Slope • Used to determine how (quickly) water flows downhill and concentrates into streams • Topographic slope can be determined from a DEM

Topographic Slope • There are three alternative sets of inputs (choose one) – Surface derivative z (dz/dx, dz/dy) – Vector with x and y components (Sx, Sy). Slope in x and y direction. – Vector with magnitude (slope) and direction (aspect) (S, )

Arc. GIS “Slope” tool • Calculates the maximum rate of change in value from that cell to its neighbors • Calculates for each cell • Represents the rate of change of elevation for each digital elevation model (DEM) cell (slope). Slope is the first derivative of a DEM • The lower the slope value, the flatter the terrain; the higher the slope value, the steeper the terrain.

Definition of X, Y, and Z in 3 D space Z axis is the direction that elevation changes (up or down) Origin is the location of the point of interest (pixel or grid cell) X axis is the direction that X has a changing value (East-West in Arc. GIS) Y axis is the direction that Y has a changing value (North -South in Arc. GIS) X, and Y are horizontal distances Z is the vertical distance The X, Y, Z axes are at right angles to one another

Definition of Slope Run is the horizontal distance calculated using X and Y Rise is the vertical distance calculated using Z (elevation) Slope ranges (-900, +900) or (-infinity %, +infinity %)

How Slope works Run is the horizontal distance calculated using X and Y Rise is the vertical distance calculated using Z (elevation)

Pythagorean theorem Used to calculate Run where a =ΔY and b = ΔX The calculated “c” is the “Run”

Arc. GIS “Slope” tool y a b c d e f g h i Calculates slope for each cell. In this illustration, it is for Cell “e” For each cell, the Slope tool calculates the maximum of the rate of change in value from that cell to each of its eight neighbors x a y b d c e g f i h x The negative sign in front of the equations is because x increases to the right (east) and y increases to the north. Now dz/dx is + if z increases with increasing x.

Arc. GIS “Slope” tool The two equations for dz/dx and dz/dy are simplified from the first equation below. The basis for that equation is illustrated in the Figure and represents an average of central finite differences over each of the three rows of cells, with the middle row counting twice as it appears in averages on each side. a y b d c e g f i h x The negative sign in front of the equations is because we are computing uphill slope

Definition of Azimuth Y axis is the direction that Y has a changing value (South to North in Arc. GIS) This is my grid cell location X axis is the direction that X has a changing value (West to East in Arc. GIS)

Definition of Azimuth Solve for α by Inverting the Tangent Function (Arc. Tan) The other way to write Arc. Tan is Tan-1 Azimuth= Convert from radians to degrees (180/π) Azimuth is the angle between North and any desired direction you want to travel

Arc. GIS Aspect – the steepest downslope direction If I pour water on the ground, which direction does it flow? Aspect is the azimuth associated with the steepest downhill slope. Therefore, we use slopes instead of distances in the tangent function. In Arc, with grid cells it is easiest to calculate Aspect using the ratio of slopes (dz/dx) and (dz/dy). = Aspect

Example for topographic slope Mesh spacing=30 m Slope/Aspect at cell e? 30 a b 80 d c 74 e 69 g f 145. 2 o 67 h 60 63 56 i 52 48 Note that this is the slope in Uphill direction (it is a positive number) Converts slope from m/m to degrees (180/π)

Example for Aspect Mesh spacing=30 m Aspect at cell e? 30 a b o -34. 8 80 74 d e 69 g 67 h 60 c 63 f 145. 2 o 56 i 52 48 One more adjustment: The above Aspect is in the direction of increasing elevation (increasing dz/dx). We need to add 180 o to this calculated aspect to get the direction of decreasing z (i. e. , the steepest downhill slope)

The Atan function is multivalued on the full circle and only unique in a range of 180 degrees. To unambiguously determine the direction from two components you really need the atan 2 function that keeps the sign on y and x components separately. For example, let y = y component of a vector x = x component of a vector atan(x/y) gives the direction of the vector as an angle (with the ratio x/y since angle here is measured from north). But x/y is the same value if y is positive and x negative, or x positive and y negative. So once you take the ratio x/y, if you get a negative number you do not know which (y or x) was negative. A way to resolve this is angle = atan(x/y) if(0 < angle < 180 and dz/dx < 0) then aspect = angle + 180 (flip the direction because dz/dx is negative) else aspect = angle endif

D 8 steepest single flow direction (Eight Direction Pour Point Model) In a gridded system, water can only flow to one of eight adjacent cells 32 64 16 8 4 128 The direction of flow is determined by the direction of steepest descent: 1 Maximum_drop = (change_in_z-value / distance) * 100 2 This is maximum percentage drop. Defined as “Hydrologic slope” in Arc. GIS ESRI Direction encoding (Arc. GIS)

Hydrologic Slope (Flow Direction Tool) Find Direction of Steepest Descent (Arc. GIS) 30 30 Slope: 80 74 63 69 67 56 60 52 48 Slope: For diagonal direction, the denominator for slope includes square root of 2

Limitation due to 8 grid directions. ? The true flow direction follows the red arrow. However, we can only choose one of the blue arrows because we have to use one of eight adjacent cells.

The D Algorithm Tarboton, D. G. , (1997), "A New Method for the Determination of Flow Directions and Contributing Areas in Grid Digital Elevation Models, " Water Resources Research, 33(2): 309 -319. ) (http: //www. engineering. usu. edu/cee/faculty/dtarb/dinf. pdf)

The D Algorithm z 3 zo z 2 z 1 If 1 does not fit within the triangle, the angle is chosen along the steepest edge or diagonal resulting in a slope and direction equivalent to D 8

D∞ Example 30 z 3 z 4 80 z 5 74 zo 69 67 z 6 z 7 60 52 284. 9 o z 2 63 56 z 1 z 8 48 14. 9 o The tool is available at http: //hydrology. usu. edu/taudem 5/documentation. html

Automating Processes using Model Builder Using a DEM tif file as input

Elevation (m) for Upper Klamath Lake Basin, OR

Elevation Contours for Wood River Valley Watershed of Upper Klamath Lake Basin

Slope (%) for Upper Klamath Lake Basin, OR (-infinity, + infinity)

Slope (Degree) for Upper Klamath Lake Basin, OR

Aspect (Degree) for Upper Klamath Lake Basin, OR

Percentage Drop (Degree) for Upper Klamath

Flow Direction Integer raster whose values range from 1 to 255 32 64 16 8 128 1 4 2

Hillshade Hypothetical illumination of a surface by determining illumination values for each cell in a raster. It does this by setting a position for a hypothetical light source and calculating the illumination values of each cell in relation to neighboring cells. Azimuth The azimuth is the angular direction of the sun, measured from north in clockwise degrees from 0 to 360. An azimuth of 90º is east. The default azimuth is 315º (NW). Altitude The altitude is the slope or angle of the illumination source above the horizon. The units are in degrees, from 0 (on the horizon) to 90 (overhead). The default is 45 degrees.

Hillshade

Map Projection Parameters

State Plane Coordinate Systems • Like UTM but customized to minimize error within States-1930 • NAD 27 – feet • NAD 83 – feet or meters 39

State Plane Coordinate Systems • • TM for North South oriented states Lambert Conformal Conic- East West orientation Depending on TM or LCC, zones baselines and meridians positioned differently 40 Baselines are placed under the zones Southern Border and measured North from there

State Plane zones State plane zones of Minnesota (1/6 of the zone width) and details of standard parallel placement (Lambert Conformal Conic)

State Plane Coordinate System-overview • Used by state/local governments – Units are feet (now appearing in meters) • Panhandle of Alaska: oblique Mercator (at angle) • N/S states (Missouri): Transverse Mercator projection • E/W states: Lambert Conformal Conic for an East-West State (California, Nebraska) • Doesn’t make sense? – States divided into smaller zones – Distortion is extremely small for local area • Coordinates: SW corner is 0, 0. • Mapping large areas (over >1 zone) is trouble!

What are the coordinates of the origin (fo, lo) and the corresponding (Xo, Yo) ? (fo, lo) = (39° 50’ N, 100° W) (lo, fo) = (100° W, 39° 50’ N) (Xo, Yo) = (1640416. 667, 0. 000) Don’t get confused. “X” is always associated with Longitude, “Y” with Latitude

Zones in the Texas State Plane

Universal Transverse Mercator Coordinate System • Line of intersection at a central meridian. • Distances measured with respect to the central meridian and from the equator in meters • 60 zones around the Earth from East to West • False Northing and False Easting 46

Universal Transverse Mercator Coordinate System • False Easting to avoid negative values: Central meridian is denoted as 500, 000 not “ 0” • Measurements are shown as positive • So if a point is 11, 254 meters west of the central meridian, it is shown as 500, 00011, 254 or 488, 746. • False northing. In the southern hemisphere, 10, 000 is the value assigned to the equator. So something 120, 000 meters south of the equator is depicted as 9, 880, 000 47

Word in UTM Grids UTM zones – lower 48 states

Universal Transverse Mercator (UTM Zone 14) Elevation (feet)

UTM (Location of N at Memorial football stadium) (Xo, Yo) = (693, 500 4, 521, 383) Meters Base map (streets) is from Arc. GIS. com

Arc. GIS. Com ready to use maps including elevation services http: //www. arcgis. com/features/maps/earth. html Elevation Land Cover Soils

Elevation Services http: //elevation. arcgis. com/arcgis/services

Summary Concepts • The elevation surface represented by a grid digital elevation model is used to derive slope important for surface flow • The eight direction pour point model approximates the surface flow using eight discrete grid directions. • The D vector surface flow model approximates the surface flow as a flow vector from each grid cell apportioned between down slope grid cells.