CEE 795 Water Resources Modeling and GIS

Lecture 4: Spatial Fields and DEM Processing (some material from Dr. David Maidment, University of Texas and Dr. David Tarboton, Utah State University) February 6, 2006 Learning Objectives:

• Demonstrate the concepts of spatial fields as a way to represent geographical information • Use raster and vector representations of spatial fields • Perform raster calculations in hydrology • Perform raster based watershed delineation from digital elevation models

Handouts: Assignments: Exercise #3

Vector and Raster Representation of Spatial Fields

Vector Raster

Numerical representation of a spatial surface ( field ) Grid TIN Contour and flowline

Six approximate representations of a field used in GIS Regularly spaced sample points Irregularly spaced sample points Rectangular Cells Irregularly shaped polygons Triangulated Irregular Network (TIN) Polylines/Contours from Longley, P. A., M. F. Goodchild, D. J. Maguire and D. W. Rind, (2001), Geographic Information Systems and Science, Wiley, 454 p.

A

grid

defines geographic space as a matrix of identically-sized square cells. Each cell holds a numeric value that measures a geographic attribute (like elevation) for that unit of space.

The grid data structure

• Grid size is defined by extent , spacing no data value information – Number of rows, number of column – Cell sizes (X and Y) – Top, left , bottom and right coordinates and • Grid values – Real (floating decimal point) – Integer (may have associated attribute table)

Definition of a Grid

Cell size Number of rows (X,Y) Number of Columns NODATA cell

Points as Cells

Line as a Sequence of Cells

Polygon as a Zone of Cells

NODATA Cells

Cell Networks

Grid Zones

Floating Point Grids

Continuous data surfaces using floating point or decimal numbers

Value attribute table for categorical (integer) grid data Attributes of grid zones

Raster Sampling

from Michael F. Goodchild. (1997) Rasters,

NCGIA Core Curriculum in GIScience

, http://www.ncgia.ucsb.edu/giscc/units/u055/u055.html, posted October 23, 1997

Raster Generalization

Largest share rule Central point rule

Raster Calculator

Cell by cell evaluation of mathematical functions

7 5 6 6 2 3 4 3 5 2 = 2 3

Example

Precipitation Losses (Evaporation, Infiltration) = Runoff

Runoff generation processes

Infiltration excess overland flow aka Horton overland flow P P q o P f f Partial area infiltration excess overland flow P P q o f P P Saturation excess overland flow P q o P q r q s

Runoff generation at a point depends on

• Rainfall intensity or amount • Antecedent conditions • Soils and vegetation • Depth to water table (topography) • Time scale of interest

These vary spatially which suggests a spatial geographic approach to runoff estimation

Modeling infiltration excess

Empirical, e.g. SCS Curve Number method Runoff from SCS Curve Number

10 9 8 2 1 0 7 6 5 4 3 0 R  ( P  0 .

2 S ) 2 P  0 .

8 S S  1000 CN  10 2 4 6

Precipitation (in)

8 10 12

Cell based discharge mapping flow accumulation of generated runoff

Radar Precipitation grid Soil and land use grid Runoff grid from raster calculator operations implementing runoff generation formula’s Accumulation of runoff within watersheds

Raster calculation – some subtleties

= +

Resampling or interpolation (and reprojection) of inputs to target extent, cell size, and projection within region defined by analysis mask Analysis mask Analysis cell size Analysis extent

Spatial Snowmelt Raster Calculation Example The grids below depict initial snow depth and average temperature over a day for an area. 150 m 100 m 40 50 55 4 6 42 47 43 2 4 42 44 41 (a) Initial snow depth (cm) (b) Temperature ( o C) One way to calculate decrease in snow depth due to melt is to use a temperature index model that uses the formula D new  D old  m  T Here D old and D new give the snow depth at the beginning and end of a time step, T gives the temperature and m is a melt factor. Assume melt factor m = 0.5 cm/ O C/day. Calculate the snow depth at the end of the day.

Snow Depth and Temperature

100 m 150 m 40 50 55 4 6 42 47 43 2 4 42 44 Initial Snow Depth (cm) 41 Temperature ( º C)

New depth calculation using Raster Calculator

[snow100m] - 0.5 * [temp150m]

The Result

38 52 41 39 • Outputs are on 150 m grid.

• How were values obtained ?

Nearest Neighbor Resampling with Cellsize Maximum of Inputs 40 50 42 47 55 43 42 44 41 40-0.5*4 = 38 55-0.5*6 = 52 42-0.5*2 = 41 41-0.5*4 = 39 4 6 38 41 52 39 2 4

Scale issues in interpretation of measurements and modeling results a) Extent

The scale triplet

b) Spacing c) Support From: Blöschl, G., (1996), Scale and Scaling in Hydrology, Habilitationsschrift, Weiner Mitteilungen Wasser Abwasser Gewasser, Wien, 346 p.

From: Blöschl, G., (1996), Scale and Scaling in Hydrology, Habilitationsschrift, Weiner Mitteilungen Wasser Abwasser Gewasser, Wien, 346 p.

Spatial analyst options for controlling the scale of the output

Extent Spacing & Support

Raster Calculator “Evaluation” of temp150

4 4 6 6

6

2 2 2 4 4 4 4 4

Nearest neighbor to the E and S has been resampled to obtain a 100 m temperature grid.

Raster calculation with options set to 100 m grid

[snow100m] - 0.5 * [temp150m]

38 47 52 41 41 45 42 41 39 • Outputs are on 100 m grid as desired.

• How were these values obtained ?

100 m cell size raster calculation 40 50 42 47 55 43 42 44 41

4 2 2 4 2 6 4 4 6 6 4 4 4

40-0.5*4 = 38 50-0.5*6 = 47 55-0.5*6 = 52 42-0.5*2 = 41 47-0.5*4 = 45 43-0.5*4 = 41 42-0.5*2 = 41 44-0.5*4 = 42 41-0.5*4 = 39 38 41 41 47 45 42 52 41 39 Nearest neighbor values resampled to 100 m grid used in raster calculation

What did we learn?

• Spatial analyst automatically uses nearest neighbor resampling • The scale (extent and cell size) can be set under options • What if we want to use some other form of interpolation?

Interpolation

Estimate values between known values.

A set of spatial analyst functions that predict values for a surface from a limited number of sample points creating a continuous raster.

Apparent improvement in resolution may not be justified

Interpolation methods • Nearest neighbor • Inverse distance weight • Bilinear interpolation z   r i 1 z i z  ( a  bx )( c  dy ) • Kriging (best linear z unbiased estimator)   w i z i • Spline z   c i x e i y e i

Nearest Neighbor “Thiessen” Polygon Interpolation Spline Interpolation

Spatial Surfaces used in Hydrology

Elevation Surface — the ground surface elevation at each point

3-D detail of the Tongue river at the WY/Mont border from LIDAR. Roberto Gutierrez University of Texas at Austin

Topographic Slope

• Defined or represented by one of the following – Surface derivative  z (dz/dx, dz/dy) – Vector with x and y components (S x , S y ) – Vector with magnitude ( slope ) and direction ( aspect ) (S,  )

Standard Slope Function a d g b e h i c f dz dx  (a  2d  g) (c  2f  i) 8 * x_mesh_spa cing dz dy  (a  2b  c) (g  2h  i) 8 * y_mesh_sp acing rise run  dz dx 2    dz dy   2 atan   dz dz / / dx dy  

Aspect – the steepest downslope direction dz dy dz dx atan   dz dz / dx / dy  

30 a

80

d

69

g

60

b

74

e

67

h

52

c

63

f 145.2

o

56

Example dz (a  2d  g) (c  2f  i)  dx 8 * x_mesh_spa cing ( 80  2 * 69  60 )  ( 63  2 * 56  24 )  8 * 30  0 .

229 i

48

dz dy  (g  2h  i) (a  2b  c) 8 * y_mesh_sp acing ( 60  2 * 52  48 )  ( 80  2 * 74  63 )  Slope  0 .

229  0 .

401 2  0 .

329 2   0 .

329 8 * 30 atan ( 0 .

401 )  21 .

8 o Aspect  atan 0 .

229  0 .

329   34 .

8 o  180 o  145 .

2 o

Hydrologic Slope - Direction of Steepest Descent 30 30

80 74 63 80 74 63 69 67 56 69 67 56 60 52 48

Slope: 67  48  0 .

45 30 2 ArcHydro Page 70

60 52 48

67  52  0 .

50 30

Eight Direction Pour Point Model

32 16 8 64 128 4 1 2 ESRI Direction encoding

ArcHydro Page 69

Limitation due to 8 grid directions.

?

Length on Meridians and Parallels (Lat, Long) = ( f , l ) Length on a Meridian: AB = R e Df (same for all latitudes) Length on a Parallel: CD = R Dl  R e Dl Cos (varies with latitude) f R e R Dl R Df R e C B A D

Example: What is the length of a 1º increment along on a meridian and on a parallel at 30N, 90W?

Radius of the earth = 6370 km.

Solution: • A 1º angle has first to be converted to radians p radians = 180 º, so 1º = p /180 = 3.1416/180 = 0.0175 radians • For the meridian, D L = R e Df  6370 * 0.0175  111 km • For the parallel,  D L = R e Dl Cos f  6370 * 0.0175 * Cos   96.5 km • Parallels converge as poles are approached 30

Summary Concepts

• Grid (raster) data structures represent surfaces as an array of grid cells • Raster calculation involves algebraic like operations on grids • Interpolation and Generalization is an inherent part of the raster data representation

Summary Concepts (2)

• The elevation surface represented by a grid digital elevation model is used to derive surfaces representing other hydrologic variables of interest such as – Slope – Drainage area (more details in later classes) – Watersheds and channel networks (more details in later classes)

Summary Concepts (3)

• The eight direction pour point model approximates the surface flow using eight discrete grid directions.