NCSU GIS 714:
Geospatial Computing and Simulations

Introduction to surface water simulation

Resources:

For animating in GRASS GIS 7, see the video instructions or instructions for the Spatio-temporal data handling and visualization in GRASS GIS workshop

For more practice in watershed analysis see Flow and watershed analysis in GIS582

Start GRASS GIS

Start GRASS - click on GRASS icon or type
grass

In startup pannel set GRASS GIS Database Directory to path to datasets, for example on MS Windows, C:\Users\myname\grassdata.
For GRASS Location select nc_spm_08_grass7 (North Carolina, State Plane, meters) and for GRASS Mapset create a new mapset (called e.g. HW2_data_simulation).
Click Start GRASS session.

Change working directory to where you will keep your external files (such as color ramps, time series layers, etc.):
Settings > GRASS working environment > Change working directory > select/create any directory
or type cd (stands for change directory) into the GUI Console and hit Enter: 

cd

Now you can use the commands from the assignment requiring external text files without the need to specify the full path to the file.

Compare flow accumulation derived by different algorithms

Set region to the high resolution study area and zoom to it. Then run least cost path flow routing with SFD and MFD options and vector based flowrouting. 
g.region rural_1m -p
r.watershed -a elevation=elev_lid792_1m threshold=5000 accumulation=accum_mfd5K drainage=draindir_5K basin=basin_mfd5K
r.watershed -as elevation=elev_lid792_1m threshold=5000 accumulation=accum_sfd5K drainage=draindir_sfd5K basin=basin_sfd5K
r.flow elevation=elev_lid792_1m flowline=flowlines flowlength=flowlg_1m flowaccumulation=flowacc_1m

Extract streams, convert to vector format and display results along with the official Wake county streams (red). The example is for the MFD results, do the same for the SFD and Dinf result and visualy compare 

r.mapcalc "streams_mfd = if(abs(accum_mfd5K) > 100, 1, null())"
r.thin streams_mfd output=streams_mfd_t
r.to.vect -s streams_mfd_t out=streams_mfd_t type=line
d.rast ortho_2001_t792_1m
d.vect streams_mfd_t
Compare the values of flowaccumulation from all your runs, including the ones from data simulation assignment, in outlet with coordinates 638845.52,220085.26. (hint, use r.what). Discuss the result.

Create a map of flooded area

Create a map of flooded area for a given water level and a seed point: 
r.lake elevation=elev_lid792_1m water_level=113.5 lake=flood1 coordinates=638728,220278
d.rast elev_lid792_1m
d.rast flood1
d.out.file floodedarea

Increase water level to 113.7 m and 114.0 m and create flooded area maps at these two levels.

Hydrology: Estimating inundation extent using HAND methodology

We will use r.stream.distance and r.lake.series are addons and we need to install them first: 
g.extension r.stream.distance
g.extension r.lake.series
We will estimate inundation extent using Height Above Nearest Drainage methodology (A.D. Nobre, 2011). We will compute HAND terrain model representing the differences in elevation between each grid cell and the elevations of the flowpath-connected downslope grid cell where the flow enters the channel. First we compute the flow accumulation, drainage and streams (with threshold value of 100000). We convert the streams to vector for better visualization. 
g.region rast=elevation
r.watershed elevation=elevation accumulation=flowacc drainage=drainage stream=streams threshold=100000
r.to.vect input=streams output=streams type=line
Now we use r.stream.distance with output parameter difference to compute new raster where each cell is the elevation difference between the cell and the the cell on the stream where the cell drains. 
r.stream.distance stream_rast=streams direction=drainage elevation=elevation method=downstream difference=above_stream
Before we compute the inundation, we will look at how r.lake works. We compute a lake from specified coordinate and water level: 
r.lake elevation=elevation water_level=90 lake=lake coordinates=637877,218475
Now instead of elevation raster we use the HAND raster to simulate 5-meter inundation and as the seed we specify the entire stream. 
r.lake elevation=above_stream water_level=5 lake=flood seed=streams
With addon r.lake.series we can create a series of inundation maps with rising water levels: 
r.lake.series elevation=above_stream start_water_level=0 end_water_level=5 water_level_step=0.5 output=inundation seed_raster=streams
r.lake.series creates a space-time dataset. We can use temporal modules to further work with the data. for example, we could further compute the volume and extent of flood water using t.rast.univar: 
t.rast.univar input=inundation separator=comma
Now we will export the results from r.lake.series so that we can display the animated series of inundation maps with rising water levels in MapBox GL. To do this we must export each raster as a png and get the bounding box of the computational region in WGS84. First download r.out.leaflet from https://github.com/ncsu-geoforall-lab/ or clone the plugin using git. 
git clone https://github.com/ncsu-geoforall-lab/grass-web-publishing.git
Next, find r.out.leaflet.py inside of the grass-web-publishing directory you just downloaded. Now in the GRASS terminal we will run the r.out.leaflet.py python script. 
mkdir images
python3 r.out.leaflet/r.out.leaflet.py raster="inundation_0.0,inundation_0.5,inundation_1.0,inundation_1.5,inundation_2.0,inundation_2.5,inundation_3.0,inundation_3.5,inundation_4.0,inundation_4.5,inundation_5.0" output="./images"
You will now see a few new files in the images directory
  • data_file.csv
  • data_file.js
  • New pngs
To add the newly output pngs to a web map continue by following the instructions found in the next tutorial here https://github.com/ncsu-geoforall-lab/grass-mapbox-tutorial