[GRASS-SVN] r62151 - in grass-addons/grass6/raster: . r.landscape.evol
svn_grass at osgeo.org
svn_grass at osgeo.org
Wed Oct 1 10:24:19 PDT 2014
Author: isaacullah
Date: 2014-10-01 10:24:19 -0700 (Wed, 01 Oct 2014)
New Revision: 62151
Added:
grass-addons/grass6/raster/r.landscape.evol/
grass-addons/grass6/raster/r.landscape.evol/585d862d.gif
grass-addons/grass6/raster/r.landscape.evol/Makefile
grass-addons/grass6/raster/r.landscape.evol/de331ef.gif
grass-addons/grass6/raster/r.landscape.evol/m100fb7e.gif
grass-addons/grass6/raster/r.landscape.evol/m11de82c.gif
grass-addons/grass6/raster/r.landscape.evol/m2c6cce6a.gif
grass-addons/grass6/raster/r.landscape.evol/m2f9c13ec.gif
grass-addons/grass6/raster/r.landscape.evol/m8e0f3ca.gif
grass-addons/grass6/raster/r.landscape.evol/r.landscape.evol
grass-addons/grass6/raster/r.landscape.evol/r.landscape.evol.html
Log:
Update to latest version of r.landscape.evol, with improved speeds and ability to read external climate file. This update offers significant improvement, better conformance to grass python coding standards, and does not break backwards compatibility with previous versions. Note that I've also decided to move r.landscape.evol out of the "LandDyn" subdirectory to make it possible to install with g.extension. This was a request from the user-base who told me it was too hard to find and install. It should be much easier now.
Added: grass-addons/grass6/raster/r.landscape.evol/585d862d.gif
===================================================================
(Binary files differ)
Property changes on: grass-addons/grass6/raster/r.landscape.evol/585d862d.gif
___________________________________________________________________
Added: svn:mime-type
+ application/octet-stream
Added: grass-addons/grass6/raster/r.landscape.evol/Makefile
===================================================================
--- grass-addons/grass6/raster/r.landscape.evol/Makefile (rev 0)
+++ grass-addons/grass6/raster/r.landscape.evol/Makefile 2014-10-01 17:24:19 UTC (rev 62151)
@@ -0,0 +1,7 @@
+MODULE_TOPDIR = ../..
+
+PGM=r.landscape.evol
+
+include $(MODULE_TOPDIR)/include/Make/Script.make
+
+default: script
Added: grass-addons/grass6/raster/r.landscape.evol/de331ef.gif
===================================================================
(Binary files differ)
Property changes on: grass-addons/grass6/raster/r.landscape.evol/de331ef.gif
___________________________________________________________________
Added: svn:mime-type
+ application/octet-stream
Added: grass-addons/grass6/raster/r.landscape.evol/m100fb7e.gif
===================================================================
(Binary files differ)
Property changes on: grass-addons/grass6/raster/r.landscape.evol/m100fb7e.gif
___________________________________________________________________
Added: svn:mime-type
+ application/octet-stream
Added: grass-addons/grass6/raster/r.landscape.evol/m11de82c.gif
===================================================================
(Binary files differ)
Property changes on: grass-addons/grass6/raster/r.landscape.evol/m11de82c.gif
___________________________________________________________________
Added: svn:mime-type
+ application/octet-stream
Added: grass-addons/grass6/raster/r.landscape.evol/m2c6cce6a.gif
===================================================================
(Binary files differ)
Property changes on: grass-addons/grass6/raster/r.landscape.evol/m2c6cce6a.gif
___________________________________________________________________
Added: svn:mime-type
+ application/octet-stream
Added: grass-addons/grass6/raster/r.landscape.evol/m2f9c13ec.gif
===================================================================
(Binary files differ)
Property changes on: grass-addons/grass6/raster/r.landscape.evol/m2f9c13ec.gif
___________________________________________________________________
Added: svn:mime-type
+ application/octet-stream
Added: grass-addons/grass6/raster/r.landscape.evol/m8e0f3ca.gif
===================================================================
(Binary files differ)
Property changes on: grass-addons/grass6/raster/r.landscape.evol/m8e0f3ca.gif
___________________________________________________________________
Added: svn:mime-type
+ application/octet-stream
Added: grass-addons/grass6/raster/r.landscape.evol/r.landscape.evol
===================================================================
--- grass-addons/grass6/raster/r.landscape.evol/r.landscape.evol (rev 0)
+++ grass-addons/grass6/raster/r.landscape.evol/r.landscape.evol 2014-10-01 17:24:19 UTC (rev 62151)
@@ -0,0 +1,642 @@
+#!/usr/bin/python
+
+############################################################################
+#
+# MODULE: r.landscape.evol.py
+# AUTHOR(S): Isaac Ullah and Michael Barton
+# COPYRIGHT: (C) 2012 GRASS Development Team/Isaac Ullah
+#
+# description: Simulates the cumulative effect of erosion and deposition on a landscape over time. This module uses appropriate flow on different landforms by default; however, singular flow regimes can be chosen by manipulating the cutoff points. This module requires GRASS 6.4 or greater. THIS SCRIPT WILL PRODUCE MANY TEMPORARY MAPS AND REQUIRES A LOT OF FREE FILE SPACE!
+
+# This program is free software; you can redistribute it and/or modify
+# it under the terms of the GNU General Public License as published by
+# the Free Software Foundation; either version 2 of the License, or
+# (at your option) any later version.
+#
+# This program is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+# GNU General Public License for more details.
+#
+#############################################################################/
+#%Module
+#% description: Simulates the cumulative effect of erosion and deposition on a landscape over time. This module uses appropriate flow on different landforms by default; however, singular flow regimes can be chosen by manipulating the cutoff points. This module requires GRASS 6.4 or greater. THIS SCRIPT WILL PRODUCE MANY TEMPORARY MAPS AND REQUIRES A LOT OF FREE FILE SPACE!
+#%End
+#%option
+#% key: elev
+#% type: string
+#% gisprompt: old,cell,raster
+#% description: Input elevation map (DEM of surface)
+#% required : yes
+#%end
+#%option
+#% key: initbdrk
+#% type: string
+#% gisprompt: old,cell,raster
+#% description: Bedrock elevations map (DEM of bedrock)
+#% answer:
+#% required : yes
+#%end
+#%option
+#% key: K
+#% type: string
+#% gisprompt: old,cell,raster
+#% description: Soil erodability index (K factor) map or constant
+#% answer: 0.42
+#% required : no
+#% guisection: Landscape Evolution
+#%end
+#%option
+#% key: sdensity
+#% type: string
+#% gisprompt: old,cell,raster
+#% description: Soil density map or constant [T/m3] for conversion from mass to volume
+#% answer: 1.2184
+#% required : no
+#% guisection: Landscape Evolution
+#%end
+#%option
+#% key: Kt
+#% type: double
+#% description: Stream transport efficiency variable (0.001 for a soft substrate, 0.0001 for a normal substrate, 0.00001 for a hard substrate, 0.000001 for a very hard substrate)
+#% answer: 0.0001
+#% options: 0.001,0.0001,0.00001,0.000001
+#% required : no
+#% guisection: Landscape Evolution
+#%end
+#%option
+#% key: loadexp
+#% type: double
+#% description: Stream transport type variable (1.5 for mainly bedload transport, 2.5 for mainly suspended load transport)
+#% answer: 1.5
+#% options: 1.5,2.5
+#% required : no
+#% guisection: Landscape Evolution
+#%end
+#%option
+#% key: kappa
+#% type: double
+#% description: Hillslope diffusion (Kappa) rate map or constant [m/kyr]
+#% answer: 1
+#% required : no
+#% guisection: Landscape Evolution
+#%end
+#%option
+#% key: C
+#% type: string
+#% gisprompt: old,cell,raster
+#% description: Landcover index (C factor) map or constant
+#% answer: 0.005
+#% required : no
+#% guisection: Landscape Evolution
+#%end
+#%option
+#% key: rain
+#% type: string
+#% gisprompt: old,cell,raster
+#% description: Precip totals for the average storm [mm] (or path to climate file of comma separated values of "rain,R,storms,stormlength", with a new line for each year of the simulation)
+#% answer: 20.61
+#% guisection: Climate
+#%end
+#%option
+#% key: R
+#% type: string
+#% description: Rainfall (R factor) constant (AVERAGE FOR WHOLE MAP AREA) (or path to climate file of comma separated values of "rain,R,storms,stormlength", with a new line for each year of the simulation)
+#% answer: 4.54
+#% guisection: Climate
+#%end
+#%option
+#% key: storms
+#% type: string
+#% description: Number of storms per year (integer) (or path to climate file of comma separated values of "rain,R,storms,stormlength", with a new line for each year of the simulation)
+#% answer: 25
+#% guisection: Climate
+#%end
+#%option
+#% key: stormlength
+#% type: string
+#% description: Average length of the storm [h] (or path to climate file of comma separated values of "rain,R,storms,stormlength", with a new line for each year of the simulation)
+#% answer: 24.0
+#% guisection: Climate
+#%end
+#%option
+#% key: speed
+#% type: double
+#% description: Average velocity of flowing water in the drainage [m/s]
+#% answer: 1.4
+#% required : no
+#% guisection: Hydrology
+#%end
+#%option
+#% key: cutoff1
+#% type: double
+#% description: Flow accumulation breakpoint value for shift from diffusion to overland flow
+#% answer: 0.65
+#% required : no
+#% guisection: Hydrology
+#%end
+#%option
+#% key: cutoff2
+#% type: double
+#% description: Flow accumulation breakpoint value for shift from overland flow to rill/gully flow (if value is the same as cutoff1, no sheetwash procesess will be modeled)
+#% answer: 2.25
+#% required : no
+#% guisection: Hydrology
+#%end
+#%option
+#% key: cutoff3
+#% type: double
+#% description: Flow accumulation breakpoint value for shift from rill/gully flow to stream flow (if value is the same as cutoff2, no rill procesess will be modeled)
+#% answer: 7
+#% required : no
+#% guisection: Hydrology
+#%end
+#%option
+#% key: smoothing
+#% type: string
+#% description: Amount of additional smoothing (answer "no" unless you notice large spikes in the erdep rate map)
+#% answer: no
+#% options: no,low,high
+#% required : yes
+#%end
+#%option
+#% key: prefx
+#% type: string
+#% description: Prefix for all output maps
+#% answer: levol_
+#% required : yes
+#%end
+#%option
+#% key: outdem
+#% type: string
+#% description: Name stem for output elevation map(s) (preceded by prefix and followed by numerical suffix if more than one iteration)
+#% answer: elevation
+#% required: yes
+#%end
+#%option
+#% key: outsoil
+#% type: string
+#% description: Name stem for the output soil depth map(s) (preceded by prefix and followed by numerical suffix if more than one iteration)
+#% answer: soildepth
+#% required: yes
+#%end
+#%option
+#% key: number
+#% type: integer
+#% description: Number of iterations (cycles) to run
+#% answer: 1
+#% required : yes
+#%end
+
+
+# #%option
+# #% key: alpha ###This may be added back in if I can find a good equation for mass movement
+# #% type: integer
+# #% description: Critical slope threshold for mass movement of sediment (in degrees above horizontal)
+# #% answer: 40
+# #% required : yes
+# #%end
+
+#%flag
+#% key: p
+#% description: -p Output a vector points map with sampled values of flow accumulation and curvatures suitable for determining cutoff values. NOTE: Overrides all other output options, and exits after completion. The output vector points map will be named "PREFIX_#_randomly_sampled_points".
+#% guisection: Optional
+#%end
+#%flag
+#% key: 1
+#% description: -1 Calculate streams as 1D difference instead of 2D divergence (forces heavy incision in streams)
+#% guisection: Optional
+#%end
+#%flag
+#% key: k
+#% description: -k Keep ALL temporary maps (overides flags -drst). This will make A LOT of maps!
+#% guisection: Optional
+#%end
+#%flag
+#% key: d
+#% description: -d Don't output yearly soil depth maps
+#% guisection: Optional
+#%end
+#%flag
+#% key: r
+#% description: -r Don't output yearly maps of the erosion/deposition rates ("ED_rate" map, in vertical meters)
+#% guisection: Optional
+#%end
+#%flag
+#% key: s
+#% description: -s Keep all slope maps
+#% guisection: Optional
+#%end
+#%flag
+#% key: t
+#% description: -t Keep yearly maps of the Transport Capacity at each cell ("Qs" maps)
+#% guisection: Optional
+#%end
+#%flag
+#% key: e
+#% description: -e Keep yearly maps of the Excess Transport Capacity (divergence) at each cell ("DeltaQs" maps)
+#% guisection: Optional
+#%end
+#%Option
+#% key: statsout
+#% type: string
+#% description: Name for the statsout text file (optional, if none provided, a default name will be used)
+#% required: no
+#% guisection: Optional
+#%end
+
+import sys
+import os
+import subprocess
+import math
+import tempfile
+grass_install_tree = os.getenv('GISBASE')
+sys.path.append(grass_install_tree + os.sep + 'etc' + os.sep + 'python')
+import grass.script as grass
+
+# Now define "main", our main block of code, here defined because of the way g.parser needs to be called with python codes for grass (see below)
+# m = last iteration number, o = iteration number, p = prefx, q = statsout, r = resolution of input elev map, s = master list of lists of climate data
+def main(m, o, p, q, r, s):
+ #get the process id to tag any temporary maps we make for easy clean up in the loop
+ pid = os.getpid()
+ # Get variables from user input
+ smoothing = os.getenv("GIS_OPT_smoothing").strip()
+ years = os.getenv("GIS_OPT_number")
+ initbdrk = os.getenv("GIS_OPT_initbdrk")
+ outdem = os.getenv("GIS_OPT_outdem")
+ outsoil = os.getenv("GIS_OPT_outsoil")
+ K = os.getenv("GIS_OPT_K")
+ sdensity = os.getenv("GIS_OPT_sdensity")
+ C = os.getenv("GIS_OPT_C")
+ kappa = os.getenv("GIS_OPT_kappa")
+ cutoff1 = os.getenv("GIS_OPT_cutoff1")
+ cutoff2 = os.getenv("GIS_OPT_cutoff2")
+ cutoff3 = os.getenv("GIS_OPT_cutoff3")
+ # these variables come in as a list of lists, so let's get this year's numbers out of them.
+ rain = s[0][m]
+ R = s[1][m]
+ storms = s[2][m]
+ stormlengthsecs = float(s[3][m])*3600.00 # number of seconds in the storm
+ stormtimet = stormlengthsecs / (float(os.getenv("GIS_OPT_speed")) * float(r)) # number of hydrologic instants in the storm
+# timet = stormlengthsecs/stormtimet # length of a single hydrologic instant in seconds, currently unused, but might be important in future versions
+ Kt = os.getenv("GIS_OPT_Kt")
+ loadexp = os.getenv("GIS_OPT_loadexp")
+ # Make some variables for temporary map names, labeled different depending on if we keep them or not
+ if ( os.getenv("GIS_FLAG_k") == "1" ):
+ aspect = '%saspect%04d' % (p, o)
+ flowacc = '%sflowacc%04d' % (p, o)
+ flacclargenums = '%sflowacc_largenums%04d' % (p, o)
+ flowdir = '%sflowdir%04d' % (p, o)
+ pc = '%spc%04d' % (p, o)
+ tc = '%stc%04d' % (p, o)
+# meancurv = '%smeancurv%04d' % (p, o) # This variable might be used if bedrock weathering is ever implemented
+# rate = '%srate%04d' % (p, o) # This variable might be used if bedrock weathering is ever implemented
+ rainexcess = "%s_rainfall_excess_map_%04d"% (p, o)
+ tempnetchange1 = '%sTEMPORARY_unsmoothed_ED_rate%04d' % (p, o)
+ tempnetchange2 = '%sTEMPORARY_smoothed_ED_rate%04d' % (p, o)
+ tmperosion = '%sTEMPORARY_erosion%04d' % (p, o)
+ tmpdep = '%sTEMPORARY_deposition%04d' % (p, o)
+ else:
+ aspect = '%saspect%04d' % (pid, o)
+ flowacc = '%sflowacc%04d' % (pid, o)
+ flacclargenums = '%sflowacc_largenums%04d' % (pid, o)
+ flowdir = '%sflowdir%04d' % (pid, o)
+ pc = '%spc%04d' % (pid, o)
+ tc = '%stc%04d' % (pid, o)
+# meancurv = '%smeancurv%04d' % (pid, o) # This variable might be used if bedrock weathering is ever implemented
+# rate = '%srate%04d' % (p, o) # This variable might be used if bedrock weathering is ever implemented
+ rainexcess = "%s_rainfall_excess_map_%04d"% (pid, o)
+ tempnetchange1 = '%sTEMPORARY_unsmoothed_ED_rate%04d' % (pid, o)
+ tempnetchange2 = '%sTEMPORARY_smoothed_ED_rate%04d' % (pid, o)
+ tmperosion = '%sTEMPORARY_erosion%04d' % (pid, o)
+ tmpdep = '%sTEMPORARY_deposition%04d' % (pid, o)
+ # Make color rules for netchange maps
+ nccolors = tempfile.NamedTemporaryFile()
+ nccolors.write('100% 0 0 100\n1 blue\n0.5 indigo\n0.01 green\n0 white\n-0.01 yellow\n-0.5 orange\n-1 red\n0% 150 0 50')
+ nccolors.flush()
+ # Make color rules for soil depth maps
+ sdcolors = tempfile.NamedTemporaryFile()
+ sdcolors.write('100% 0:249:47\n20% 78:151:211\n6% 194:84:171\n0% 227:174:217')
+ sdcolors.flush()
+ # If first iteration, use input maps. Otherwise, use maps generated from previous iterations
+ if ( o == 1 ):
+ old_dem = '%s' % os.getenv("GIS_OPT_elev")
+ old_bdrk = '%s' % os.getenv("GIS_OPT_initbdrk")
+ old_soil = "%s%s_init" % (prefx, os.getenv("GIS_OPT_outsoil"))
+ grass.mapcalc('${old_soil}=${old_dem}-${old_bdrk}', overwrite = "True", quiet = "True", old_soil = old_soil, old_dem = old_dem, old_bdrk = old_bdrk)
+ else :
+ old_dem = '%s%s%04d' % (p, os.getenv("GIS_OPT_outdem"), m)
+ old_bdrk = '%s%s%04d' % (p, os.getenv("GIS_OPT_outbdrk"), m)
+ old_soil = '%s%s%04d' % (p, os.getenv("GIS_OPT_outsoil"), m)
+ #Checking for special condition of there being only one run, and setting variables accordingly (one year runs have no numbers suffixed to the output map names)
+ if ( years == '1' ):
+ slope = '%sslope' % p
+ netchange = '%sED_rate' % p
+ new_dem ='%s%s' % (p, outdem)
+ new_soil = '%s%s' % (p, outsoil)
+ else:
+ slope = '%sslope%04d' % (p, o)
+ netchange = '%sED_rate%04d' % (p, o)
+ new_dem = '%s%s%04d' % (p, outdem, o)
+ new_soil = '%s%s%04d' % (p, outsoil, o)
+ #Check to see if we are going to only output diagnostics for determing cutoff values, and act accordingly
+ if ( os.getenv("GIS_FLAG_p") == "1" ):
+ grass.message('GATHERING STATISTICS FOR DETERMINING CUTOFF VALUES\n-------------------------------------------------\n1) Calculating slope and curvatures')
+ grass.run_command('r.slope.aspect', quiet = "True", elevation = old_dem, slope = slope, pcurv = pc, tcurv = tc)
+ else:
+ grass.message('\n##################################################\n\n*************************\n Iteration %s -- ' % o + 'step 1: calculating slope\n*************************\n')
+ grass.run_command('r.slope.aspect', quiet = "True", elevation = old_dem, aspect = aspect, slope = slope)
+ if ( os.getenv("GIS_FLAG_p") == "1" ):
+ grass.message('2) Calculating map of rainfall excess')
+ else:
+ grass.message('\n*************************\n Iteration %s -- ' % o + 'step 2: calculating accumulated flow depths\n*************************\n')
+ grass.message('Calculating runoff excess rates (scaled uplsope accumulated cells)')
+
+ #to calculate rainfall excess, we are making a linear regression of C factor and the percentage of water that will leave the cell. A c-factor of 0.005 (mature woodland) will only allow 10% of the water to exit the cell, whereas a c-factor of 0.1 (bareland) will allow 98% of the water to leave the cell. Note that we multiply this by 100 because r.watershed will only allow integer amounts as input in it's 'flow' variable, and we want to maintain the accuracy. The large number flow accumulation will be divided by 100 after it is made, which brings the values back down to what they should be.
+ grass.mapcalc('${rainexcess}=int(100 * ((9.26316 * ${C}) + 0.05368))', quiet = "True", rainexcess = rainexcess, C = C)
+ if ( os.getenv("GIS_FLAG_p") == "1" ):
+ grass.message('3) Calculating accumulated flow (in numbers of upslope cells, scaled by runoff contribution')
+ #r.watershed -a elevation=hvdf flow=fad accumulation=sadasd
+ try:
+ grass.run_command('r.watershed', quiet = "True", flags = 'fa', elevation = old_dem, flow = rainexcess, accumulation = flacclargenums, drainage = flowdir, convergence = '5')
+ except:
+ grass.message("using grass 7 syntax for r.watershed instead")
+ grass.run_command('r.watershed', quiet = "True", flags = 'a', elevation = old_dem, flow = rainexcess, accumulation = flacclargenums, drainage = flowdir, convergence = '5')
+ grass.mapcalc('${flowacc}=${flacclargenums}/100', quiet = "True", flowacc = flowacc, flacclargenums = flacclargenums)
+ #again, do something different if we are only making an evaluation of cutoffs
+ if ( os.getenv("GIS_FLAG_p") == "1" ):
+ grass.message('4) Determining number of sampling points using formula: "ln(#cells_in_input_map)*100"')
+ flaccstats = grass.parse_command('r.univar', flags = 'g', map = flowacc)
+ numpts = int(math.log(int(flaccstats['n']))*100)
+ grass.message('5) Creating random points and sampling values of flow accumulation, curvatures, and slope.')
+ vout = '%s%s_randomly_sampled_points' % (p, numpts)
+ grass.run_command('r.random', quiet = "True", input = flowacc, cover = pc, n = numpts, vector_output = vout)
+ grass.run_command('v.db.renamecol', quiet = "True", map = vout, column = 'value,Flow_acc')
+ grass.run_command('v.db.renamecol', quiet = "True", map = vout, column = 'covervalue,Princ_curv')
+ grass.run_command('v.db.addcol', quiet = "True", map = vout, columns = 'Tang_curv double precision, Slope double precision')
+ grass.run_command('v.what.rast', quiet = "True", vector = vout, raster = tc, column = "Tang_curv")
+ grass.run_command('v.what.rast', quiet = "True", vector = vout, raster = slope, column = "Slope")
+ if ( os.getenv("GIS_FLAG_k") == "1" ):
+ grass.message('--Keeping the created maps (Flow Accumulation, Slope, Principle Curvature, Tangential Curvature)')
+ else:
+ grass.message('6) Cleaning up...')
+ grass.run_command('g.remove', quiet = "True", rast = slope + "," + pc + "," + tc + "," + flowacc)
+ grass.message('FINISHED. \nRandom sample points map "%s" created successfully.\n' % vout)
+ sys.exit(0)
+ grass.message('\n*************************\n Iteration %s -- ' % o + 'step 3: calculating sediment transport rates (units variable depending upon process) \n*************************\n')
+ # This step calculates the force of the flowing water at every cell on the landscape using the proper transport process law for the specific point in the flow regime. For upper hillslopes (below cutoff point 1) this done by multiplying the diffusion coeficient by the accumulated flow/cell res width. For midslopes (between cutoff 1 and 2) this is done by multiplying slope by accumulated flow with the m and n exponents set to 1. For channel catchment heads (between cutoff 2 and 3), this is done by multiplying slope by accumulated flow with the m and n exponents set to 1.6 and 1.3 respectively. For Channelized flow in streams (above cutoff 3), this is done by calculating the reach average shear stress (hydraulic radius [here estimated for a cellular landscape simply as the depth of flow] times slope times accumulated flow [cells] times gravitatiopnal acceleration of water [9806.65 newtons], all raised to the appropriate exponant for the type of transport (bedload or suspe
nded load), and then divided by the resolution. Depth of flow is calculated as a mean "instantaneous depth" during any given rain event, here estimated by the maximum depth of an idealized unit hydrograph with base equal to the duration of the storm, and area equal to the total accumulated excess rainfall during the storm. Then finally calculates the stream power or sediment carrying capacity (qs) of the water flowing at each part of the map by multiplying the reach average shear stress (channelized flow in streams) or the estimated flow force (overland flow) by the transport coeficient (estimated by R*K*C for hillslopes or kt for streams). This is a "Transport Limited" equation, however, we add some constraints on detachment by checking to see if the sediment supply has been exhausted: if the current soil depth is 0 or negative (checking for a negative value is kind of an error trap) then we make the transport coefficient small (0.000001) to simulate erosion on bedrock. Bec
ause diffusion and USPED require 2D divergence later on, we calculate these as vectors in the X and Y directions. Stream flow only needs 1D difference, so it's calulated in the direction of flow.
+ if ( os.getenv("GIS_FLAG_1") == "1" ):
+ #This is the version with 1D streams
+ qs1 = '%sQs_1D_streams%04d' % (p, o)
+ grass.mapcalc('${qs1}=(${Kt} * exp(9806.65*(((${rain}/1000)*${flowacc})/(0.595*${stormtimet}))*sin(${slope}), ${loadexp}) )', quiet = "True", qs1 = qs1, flowacc = flowacc, stormtimet = stormtimet, slope = slope, loadexp = loadexp, Kt = Kt, sdensity = sdensity)
+ qsx = "%sQsx_%04d" % (p,o)
+ qsy = "%sQsy_%04d" % (p,o)
+ grass.mapcalc("${qsx}=eval(a=(${kappa} * sin(${slope}) * cos(${aspect})), b=((${R}*${K}*${C}*${flowacc}*${res}*sin(${slope})) * cos(${aspect})), c=( (${R}*${K}*${C}*exp((${flowacc}*${res}),1.6000000)*exp(sin(${slope}),1.3000000)) * cos(${aspect})), if(${flowacc} >= ${cutoff3}, a, if(${flowacc} >= ${cutoff2} && ${flowacc} < ${cutoff3}, b, c)) )", quiet = "True", qsx = qsx, kappa = kappa, slope = slope, aspect = aspect, R = R, K = K, C =C, res = r, flowacc = flowacc, Kt = Kt, rain = rain, stormtimet = stormtimet, loadexp = loadexp, cutoff1 = cutoff1, cutoff2 = cutoff2, cutoff3 = cutoff3)
+ grass.mapcalc("${qsy}=eval(a=(${kappa} * sin(${slope}) * sin(${aspect})), b=((${R}*${K}*${C}*${flowacc}*${res}*sin(${slope})) * sin(${aspect})), c=( (${R}*${K}*${C}*exp((${flowacc}*${res}),1.6000000)*exp(sin(${slope}),1.3000000)) * sin(${aspect})), d=(${Kt} * exp(9806.65*(((${rain}/1000)*${flowacc})/(0.595*${stormtimet}))*sin(${slope}), ${loadexp}) ) * sin(${aspect}), if(${flowacc} >= ${cutoff3}, a, if(${flowacc} >= ${cutoff2} && ${flowacc} < ${cutoff3}, b, c)) )", quiet = "True", qsy = qsy, kappa = kappa, slope = slope, aspect = aspect, R = R, K = K, C =C, res = r, flowacc = flowacc, Kt = Kt, rain = rain, stormtimet = stormtimet, loadexp = loadexp, cutoff1 = cutoff1, cutoff2 = cutoff2, cutoff3 = cutoff3)
+ else:
+ #This is the normal version (with 2D streams)
+ qsx = "%sQsx_%04d" % (p,o)
+ qsy = "%sQsy_%04d" % (p,o)
+ grass.mapcalc("${qsx}=eval(a=(${kappa} * sin(${slope}) * cos(${aspect})), b=((${R}*${K}*${C}*${flowacc}*${res}*sin(${slope})) * cos(${aspect})), c=( (${R}*${K}*${C}*exp((${flowacc}*${res}),1.6000000)*exp(sin(${slope}),1.3000000)) * cos(${aspect})), d=(${Kt} * exp(9806.65*(((${rain}/1000)*${flowacc})/(0.595*${stormtimet}))*sin(${slope}), ${loadexp}) ) * cos(${aspect}), if(${flowacc} >= ${cutoff3}, a, if(${flowacc} >= ${cutoff2} && ${flowacc} < ${cutoff3}, b, if(${flowacc} >= ${cutoff1} && ${flowacc} < ${cutoff2}, c, d))) )", quiet = "True", qsx = qsx, kappa = kappa, slope = slope, aspect = aspect, R = R, K = K, C =C, res = r, flowacc = flowacc, Kt = Kt, rain = rain, stormtimet = stormtimet, loadexp = loadexp, cutoff1 = cutoff1, cutoff2 = cutoff2, cutoff3 = cutoff3)
+ grass.mapcalc("${qsy}=eval(a=(${kappa} * sin(${slope}) * sin(${aspect})), b=((${R}*${K}*${C}*${flowacc}*${res}*sin(${slope})) * sin(${aspect})), c=( (${R}*${K}*${C}*exp((${flowacc}*${res}),1.6000000)*exp(sin(${slope}),1.3000000)) * sin(${aspect})), d=(${Kt} * exp(9806.65*(((${rain}/1000)*${flowacc})/(0.595*${stormtimet}))*sin(${slope}), ${loadexp}) ) * sin(${aspect}), if(${flowacc} >= ${cutoff3}, a, if(${flowacc} >= ${cutoff2} && ${flowacc} < ${cutoff3}, b, if(${flowacc} >= ${cutoff1} && ${flowacc} < ${cutoff2}, c, d))) )", quiet = "True", qsy = qsy, kappa = kappa, slope = slope, aspect = aspect, R = R, K = K, C =C, res = r, flowacc = flowacc, Kt = Kt, rain = rain, stormtimet = stormtimet, loadexp = loadexp, cutoff1 = cutoff1, cutoff2 = cutoff2, cutoff3 = cutoff3)
+
+ grass.message('\n*************************\n Iteration %s -- ' % o + 'step 4: calculating divergence/difference of sediment transport for each process and the actual amount of erosion or deposition in vertical meters/cell/year\n*************************\n\n')
+ #Here is where we figure out the change in transport capacity, and thus the actual amount of erosion an deposition that would occur. There are two ways of doing this. On planar and convex surfaces (i.e., ridgetops, flats, hillslopes), it is better to take the 2D divergence of sediment flux (we use r.slope.aspect to calculate this), but on highly convex surfaces (i.e., in channels) it is better to take the 1D difference between one cell, and the cell that is immediately downstream from it. This all assumes that the system is always operating at Transport Capacity, or if it is not, then is still behaves as if it were (ie., that the actual differences in transported sediment between the cells would be proportional to the system operating at capacity). Thus, under this assumption, the divergence of capacity is equals to actual amount of sediment eroded/deposited.
+ #This is the way we implemnt this: First calculate, we calculate the divergence/differnce for EACH of the different flow processes on the ENTIRE map (i.e., make one map per process, difference for streams, divergence for USPED and diffusion). Then, we cut out the pieces of each of these maps that correspond to the correct landforms from each specific process (based on the user-input cutoffs in flow accumulation), and patch them together into a single map (NOTE: see output unit conversions section below to see how we get all the units to line up during this process). This counters the "boundary effect" that happens when running the differential equations for divergence across the boundary of two different flow processes. Then we may still have to run a median smoother on the patched map to get rid of any latent spikes.
+ if ( os.getenv("GIS_FLAG_1") == "1" ):
+ #This is the version with 1D streams
+ qsd1 = '%sDelta_Qs_1D_streams%04d' % (p, o)
+ grass.mapcalc('${qsd1}=if(${flowdir} == 7, (${qs1}[-1,-1]-${qs1}), if (${flowdir} == 6, (${qs1}[-1,0]-${qs1}), if (${flowdir} == 5, (${qs1}[-1,1]-${qs1}), if (${flowdir} == 4, (${qs1}[0,1]-${qs1}), if (${flowdir} == 3, (${qs1}[1,1]-${qs1}), if (${flowdir} == 2, (${qs1}[1,0]-${qs1}), if (${flowdir} == 1, (${qs1}[1,-1]-${qs1}), if (${flowdir} == 8, (${qs1}[0,-1]-${qs1}), ${qs1}))))))))', quiet = "True", qsd1 = qsd1, flowdir = flowdir, qs1 = qs1)
+ qsxdx = '%sDelta_Qsx_%04d' % (p, o)
+ qsydy = '%sDelta_Qsy_%04d' % (p, o)
+ grass.run_command('r.slope.aspect', quiet = "True", elevation = qsx, dx = qsxdx)
+ grass.run_command('r.slope.aspect', quiet = "True", elevation = qsy, dy = qsydy)
+ else:
+ #This is the normal version (with 2D streams)
+ qsxdx = '%sDelta_Qsx_%04d' % (p, o)
+ qsydy = '%sDelta_Qsy_%04d' % (p, o)
+ grass.run_command('r.slope.aspect', quiet = "True", elevation = qsx, dx = qsxdx)
+ grass.run_command('r.slope.aspect', quiet = "True", elevation = qsy, dy = qsydy)
+
+ #This is the smoothing routine. First we calculate the rate of Erosion and Deposition by converting the Delta QS of the different processes to vertical meters by dividing by the soil denisity (with apropriate constants to get into the correct units, see UNIT CONVERSION note below), and for streams, also expand from the storm to the year level. All units of this initial (temporary) ED_rate map will be in m/cell/year.
+ #OUTPUT UNIT CONVERSIONS: In the case of the diffusion equation, the output units are in verticle meters of sediment per cell per year, so these will be left alone. In the case of stream flow, the output units are kg/m2/storm, so need to multiply by 1000 to get T/m2/storm, and then divide by the soil density (T/m3) to get verticle meters of sediment/cell/storm, and will be multiplied by the number of storms/year in order to get vertical meters of sediment/cell/year. In the case of USPED, the output is in T/Ha/year, so first multiply by 0.1 to get T/m2/year and then divide by soil density (T/m3) to get verticle meters of sediment/cell/year.
+ if ( os.getenv("GIS_FLAG_1") == "1" ):
+ #This is the version with 1D streams
+ grass.mapcalc('${tempnetchange1}=if(${flowacc} >= ${cutoff3}, ((${qsd1})/(${sdensity}*1000))*(${storms}*0.25*${stormtimet}), if(${flowacc} >= ${cutoff1} && ${flowacc} < ${cutoff3}, ((${qsxdx}+${qsydy})*0.1)/${sdensity}, ${qsxdx}+${qsydy}))', quiet = "True", tempnetchange1 = tempnetchange1, qsd1 = qsd1, qsxdx = qsxdx, qsydy = qsydy, flowacc = flowacc, cutoff1 = cutoff1, cutoff3 = cutoff3, sdensity = sdensity, storms = storms, stormtimet = stormtimet)
+ else:
+ #This is the normal version (with 2D streams)
+ grass.mapcalc('${tempnetchange1}=if(${flowacc} >= ${cutoff3}, ((${qsxdx} + ${qsydy})/(${sdensity}*1000))*(${storms}*0.25*${stormtimet}), if(${flowacc} >= ${cutoff1} && ${flowacc} < ${cutoff3}, ((${qsxdx}+${qsydy})*0.1)/${sdensity}, ${qsxdx}+${qsydy}))', quiet = "True", tempnetchange1 = tempnetchange1, qsxdx = qsxdx, qsydy = qsydy, flowacc = flowacc, cutoff1 = cutoff1, cutoff3 = cutoff3, sdensity = sdensity, storms = storms, stormtimet = stormtimet)
+
+ #Make some temp maps of just erosion rate and just deposition rate so we can grab some stats from them for the soft-knee limiting filter
+ grass.message('Running soft-knee smoothing filter...')
+ grass.mapcalc('${tmperosion}=if(${tempnetchange1} < -0, ${tempnetchange1}, null())', quiet = "True", tmperosion = tmperosion, tempnetchange1 = tempnetchange1)
+ grass.mapcalc('${tmpdep}=if(${tempnetchange1} > 0, ${tempnetchange1}, null())', quiet = "True", tmpdep = tmpdep, tempnetchange1 = tempnetchange1)
+ #Grab the stats from these temp files and save them to dictionaries
+ erosstats = grass.parse_command('r.univar', flags = 'ge', percentile = '1', map = tmperosion)
+ depostats = grass.parse_command('r.univar', flags = 'ge', percentile = '99', map = tmpdep)
+ maximum = depostats['max']
+ minimum = erosstats['min']
+ erosbreak = float(erosstats['first_quartile'])
+ deposbreak = float(depostats['third_quartile'])
+ scalemin = float(erosstats['percentile_1'])
+ scalemax = float(depostats['percentile_99'])
+ #Use the stats we gathered to do some smoothing with a hi-cut and lo-cut filter (with soft-knee limiting) of the unsmoothed ED_rate map. Values from the 1st quartile of erosion to the minimum (i.e., the very large negative numbers) will be rescaled linearly from the 1st quartile to the 1st percentile value, and values from the 3rd quartile of deposition to the maximum (i.e., the very large positiive numbers) will be rescaled linearly from the 3rd quartile to the 99th percentile value. This brings any values that were really unreasonnable as originally calculated (spikes) into the range of what the maximum values should be on a normally distrubuted dataset, but does so with out a "brick wall" style of limiting, which would make all values above some cutoff equal to a theoretical maximum. By setting both maximum cutoff point AND a "soft" scaling point, this "soft-knee" style of limiting sill retains some of the original scaling at the high ends, which allows for the smooth
ed value of very high cells to still be relatively higher than values in other cells that were also above the scaling cutoff, but were not originally as high as those very high cells.
+ grass.mapcalc('${tempnetchange2}=graph(${tempnetchange1}, ${minimum},${scalemin}, ${erosbreak},${erosbreak}, ${deposbreak},${deposbreak}, ${maximum},${scalemax})', quiet = "True", tempnetchange2 = tempnetchange2, tempnetchange1 =tempnetchange1, minimum = minimum, scalemin = scalemin, erosbreak = erosbreak, deposbreak = deposbreak, maximum = maximum, scalemax = scalemax)
+ #Check if additional smoothing is requested.
+ if smoothing == "no":
+ grass.message('No additional modal smoothing was requested...')
+ grass.run_command('g.rename', quiet = "True", rast = tempnetchange2 + ',' + netchange)
+ elif smoothing == "low":
+ grass.message('Enacting additional "low" smoothing: one pass of a 3x3 modal smoothing window.')
+ grass.run_command('r.neighbors', quiet = "True", input = tempnetchange2, output = netchange, method = 'mode', size = '3')
+ elif smoothing == "high":
+ grass.message('Enacting additional "high" smoothing: one pass of a 5x5 modal smoothing window.')
+ grass.run_command('r.neighbors', quiet = "True", input = tempnetchange2, output = netchange, method = 'mode', size = '5')
+ else:
+ grass.message('There was a problem reading the median-smoothing variable, so maps will not be median-smoothed.')
+ grass.run_command('g.rename', quiet = "True", rast = tempnetchange2 + ',' + netchange)
+ #Set the netchange map colors to the rules we've provided above
+ grass.run_command('r.colors', quiet = "True", map = netchange, rules = nccolors.name)
+ #Grab the stats from these new smoothed netchange maps and save them to dictionaries (Note that the temporary erosion and deposition maps made in this step are overwriting the two temporary maps made for gathering the stats for the soft-knee limiting filter)
+ grass.mapcalc('${tmperosion}=if(${netchange} < -0, ${netchange}, null())', quiet = "True", overwrite = "True", tmperosion = tmperosion, netchange = netchange)
+ grass.mapcalc('${tmpdep}=if(${netchange} > 0, ${netchange}, null())', quiet = "True", overwrite = "True", tmpdep = tmpdep, netchange = netchange)
+ erosstats1 = grass.parse_command('r.univar', flags = 'ge', map = tmperosion)
+ depostats1 = grass.parse_command('r.univar', flags = 'ge', map = tmpdep)
+
+ grass.message('\n*************************\n Iteration %s -- ' % o + 'step 5: calculating terrain evolution and new soil depths\n *************************\n\n')
+ #Set up a temp dem, and then do initial addition of ED change to old DEM. This mapcalc statement first checks the amount of erodable soil in a given cell against the amount of erosion calculated, and keeps the cell from eroding past this amount (if there is soil, then if the amount of erosion is more than the amount of soil, just remove all the soil and stop, else remove the amount of caclulated erosion. It also runs an error catch that checks to make sure that soil depth is not negative (could happen, I suppose), and if it is, corrects it). Finally, do patch-job to catch the shrinking edge problem (the edge cells have no upstream cell, so get turned null in the calculations in step 4)
+ grass.mapcalc('${new_dem}=eval(x=if(${old_soil} > 0.0 && (-1*${netchange}) <= ${old_soil}, ${netchange}, if((-1*${netchange}) > ${old_soil}, (-1*${old_soil}), 0)), y=(${old_dem} + x), if(isnull(y), ${old_dem}, y))', quiet = "True", new_dem = new_dem, old_soil = old_soil, old_dem = old_dem, netchange = netchange)
+ #Set colors for elevation map to match other dems
+ grass.run_command('r.colors', quiet = "True", map = new_dem, rast = os.getenv("GIS_OPT_elev"))
+ grass.mapcalc('${new_soil}=if ((${new_dem} - ${initbdrk}) < 0, 0, (${new_dem} - ${initbdrk}))', quiet = "True", new_soil = new_soil, new_dem = new_dem, initbdrk = initbdrk)
+ grass.run_command('r.colors', quiet = "True", map = new_soil, rules = sdcolors.name)
+ grass.message('\n*************************\n Iteration %s -- ' % o + 'step 6: writing stats to output file\n *************************\n\n')
+ #Finish gathering stats (just need the soil depth stats now)
+ soilstats = grass.parse_command('r.univar', flags = 'ge', map = new_soil, percentile = '99')
+ #Write stats to a new line in the stats file
+ #HEADER of the file should be: ',,Mean Values,,,,Standard Deviations,,,,Totals,,,Additional Stats\nIteration,,Mean Erosion,Mean Deposition,Mean Soil Depth,,Standard Deviation Erosion,Standard Deviation Deposition,Standard Deviation Soil Depth,,Total Sediment Eroded,Total Sediment Deposited,,Minimum Erosion,First Quartile Erosion,Median Erosion,Third Quartile Erosion,Maximum Erosion,Original Un-smoothed Maximum Erosion,,Minimum Deposition,First Quartile Deposition,Median Deposition,Third Quartile Deposition,Maximum Deposition,Original Un-smoothed Maximum Deposition,,Minimum Soil Depth,First Quartile Soil Depth,Median Soil Depth,Third Quartile Soil Depth,Maximum Soil Depth'
+ grass.message('Outputing stats to textfile: ' + q)
+ f.write('\n%s' % o + ',,' + erosstats1['mean'] + ',' + depostats1['mean'] + ',' + soilstats['mean'] + ',,' + erosstats1['stddev'] + ',' + depostats1['stddev'] + ',' + soilstats['stddev'] + ',,' + erosstats1['sum'] + ',' + depostats1['sum'] + ',,' + erosstats1['max'] + ',' + erosstats1['third_quartile'] + ',' + erosstats1['median'] + ',' + erosstats1['first_quartile'] + ',' + erosstats1['min'] + ',' + minimum + ',,' + depostats1['min'] + ',' + depostats1['first_quartile'] + ',' + depostats1['median'] + ',' + depostats1['third_quartile'] + ',' + depostats1['max'] + ',' + maximum + ',,' + soilstats['min'] + ',' + soilstats['first_quartile'] + ',' + soilstats['median'] + ',' + soilstats['third_quartile'] + ',' + soilstats['max'])
+
+ #Clean up temporary maps
+ if os.getenv("GIS_FLAG_k") == "1":
+ grass.message('\nTemporary maps will NOT be deleted!!!!\n')
+ else:
+ grass.message('\nCleaning up temporary maps...\n\n')
+ #first remove all the easy temporary maps labeled with "pid"
+ try:
+ grass.run_command("g.mremove", quiet = "True", flags = 'f', rast = "%s*" % pid)
+ except:
+ grass.run_command("g.mremove", quiet = "True", flags = 'f', type = 'rast', pattern = '%s*' % pid)
+ #now check all the flag options, and build a list of maps to delete
+ mapstoremove = []
+ if os.getenv("GIS_FLAG_s") == "1":
+ grass.message('Keeping Slope map.')
+ else:
+ mapstoremove.append(slope)
+ if os.getenv("GIS_FLAG_d") == "1":
+ grass.message('Not keeping Soil Depth map.')
+ mapstoremove.append(old_soil)
+ #check if this is the last year and remove the "new-soil" map too
+ if ( o == int(os.getenv("GIS_OPT_number"))):
+ mapstoremove.append(new_soil)
+ else:
+ #check if this is the first year, and if so, remove the temporary "soildepths_init" map
+ if ( o <= 1 ):
+ mapstoremove.append("%s%s_init" % (prefx, os.getenv("GIS_OPT_outsoil")))
+ if ( os.getenv("GIS_FLAG_e") == "1" ):
+ grass.message('Keeping Excess Transport Capacity (divergence) maps for all processes.')
+ else:
+ mapstoremove.extend([qsxdx, qsydy])
+ if ( os.getenv("GIS_FLAG_1") == "1" ):
+ mapstoremove.append(qsd1)
+ if ( os.getenv("GIS_FLAG_t") == "1" ):
+ grass.message('Keeping Transport Capacity maps for all processes.')
+ else:
+ mapstoremove.extend([qsx, qsy])
+ if ( os.getenv("GIS_FLAG_1") == "1" ):
+ mapstoremove.append(qs1)
+ if ( os.getenv("GIS_FLAG_r") == "1" ):
+ grass.message('Not keeping an Erosion and Deposition rate map.')
+ mapstoremove.append(netchange)
+ if len(mapstoremove) == 0:
+ pass
+ else:
+ grass.run_command('g.remove', quiet = "True", rast = ','.join(mapstoremove))
+ sdcolors.close()
+ nccolors.close()
+ grass.message('\n*************************\nDone with Iteration %s ' % o + '\n*************************\n')
+
+#Here is where the code in "main" actually gets executed. This way of programming is neccessary for the way g.parser needs to run.
+if __name__ == "__main__":
+ if ( len(sys.argv) <= 1 or sys.argv[1] != "@ARGS_PARSED@" ):
+ os.execvp("g.parser", [sys.argv[0]] + sys.argv)
+ else:
+ # Set up some basic variables
+ years = os.getenv("GIS_OPT_number")
+ prefx = os.getenv("GIS_OPT_prefx")
+ #these values could be read in from a climate file, so check that, and act accordingly. Either way, the result will be some lists with the same number of entries as there are iterations.
+ rain2 = []
+ try:
+ rain1 = float(os.getenv("GIS_OPT_rain"))
+ for year in range(int(years)):
+ rain2.append(rain1)
+ except:
+ with open(os.getenv("GIS_OPT_rain"), 'rU') as f:
+ for line in f:
+ rain2.append(line.split(",")[0])
+ #check for text header and remove if present
+ try:
+ float(rain2[0])
+ except:
+ del rain2[0]
+ #throw a warning if there aren't enough values in the column
+ if len(rain2) != int(years):
+ grass.fatal("Number of rows of rainfall data in your climate file\n do not match the number of iterations you wish to run.\n Please ensure that these numbers match and try again")
+ sys.exit(1)
+ R2 = []
+ try:
+ R1 = float(os.getenv("GIS_OPT_R"))
+ for year in range(int(years)):
+ R2.append(R1)
+ except:
+ with open(os.getenv("GIS_OPT_R"), 'rU') as f:
+ for line in f:
+ R2.append(line.split(",")[1])
+ #check for text header and remove if present
+ try:
+ float(R2[0])
+ except:
+ del R2[0]
+ #throw a warning if there aren't enough values in the column
+ if len(R2) != int(years):
+ grass.fatal("Number of rows of R-Factor data in your climate file\n do not match the number of iterations you wish to run.\n Please ensure that these numbers match and try again")
+ sys.exit(1)
+ storms2 = []
+ try:
+ storms1 = float(os.getenv("GIS_OPT_storms"))
+ for year in range(int(years)):
+ storms2.append(storms1)
+ except:
+ with open(os.getenv("GIS_OPT_storms"), 'rU') as f:
+ for line in f:
+ storms2.append(line.split(",")[2])
+ #check for text header and remove if present
+ try:
+ float(storms2[0])
+ except:
+ del storms2[0]
+ #throw a warning if there aren't enough values in the column
+ if len(storms2) != int(years):
+ grass.fatal("Number of rows of storm frequency data in your climate file\n do not match the number of iterations you wish to run.\n Please ensure that these numbers match and try again")
+ sys.exit(1)
+ stormlength2 = []
+ try:
+ stormlength1 = float(os.getenv("GIS_OPT_stormlength"))
+ for year in range(int(years)):
+ stormlength2.append(stormlength1)
+ except:
+ with open(os.getenv("GIS_OPT_stormlength"), 'rU') as f:
+ for line in f:
+ stormlength2.append(line.split(",")[3])
+ #check for text header and remove if present
+ try:
+ float(stormlength2[0])
+ except:
+ del stormlength2[0]
+ #throw a warning if there aren't enough values in the column
+ if len(stormlength2) != int(years):
+ grass.fatal("Number of rows of storm length data in your climate file\n do not match the number of iterations you wish to run.\n Please ensure that these numbers match and try again")
+ sys.exit(1)
+ #Now gather these four lists into one master list, to make it easier to pass on to main()
+ masterlist = [rain2,R2,storms2,stormlength2]
+ #Make the statsout file with correct column headers
+ if os.getenv("GIS_OPT_statsout") == "":
+ env = grass.gisenv()
+ mapset = env['MAPSET']
+ statsout = '%s_%slsevol_stats.csv' % (mapset, prefx)
+ else:
+ statsout = os.getenv("GIS_OPT_statsout")
+ if os.path.isfile(statsout):
+ f = file(statsout, 'a')
+ else:
+ f = file(statsout, 'wt')
+ f.write('These statistics are in units of vertical meters (depth) per cell\n,,Mean Values,,,,Standard Deviations,,,,Totals,,,Additional Stats\nIteration,,Mean Erosion,Mean Deposition,Mean Soil Depth,,Standard Deviation Erosion,Standard Deviation Deposition,Standard Deviation Soil Depth,,Total Sediment Eroded,Total Sediment Deposited,,Minimum Erosion,First Quartile Erosion,Median Erosion,Third Quartile Erosion,Maximum Erosion,Original Un-smoothed Maximum Erosion,,Minimum Deposition,First Quartile Deposition,Median Deposition,Third Quartile Deposition,Maximum Deposition,Original Un-smoothed Maximum Deposition,,Minimum Soil Depth,First Quartile Soil Depth,Median Soil Depth,Third Quartile Soil Depth,Maximum Soil Depth')
+ if ( os.getenv("GIS_FLAG_p") == "1" ):
+ grass.message('Making sample points map for determining cutoffs.')
+ else:
+ grass.message('\n##################################################\n##################################################\n\n STARTING SIMULATION\n\nBeginning iteration sequence. This may take some time.\nProcess is not finished until you see the message: \'Done with everything\'\n _____________________________________________________________\n_____________________________________________________________\n')
+ grass.message("Total number of iterations to be run is %s" % years)
+ #Get the region settings
+ region1 = grass.region()
+ # This is the loop!
+ for x in range(int(years)):
+ grass.message("Iteration = %s" % (x + 1))
+ main(x, (x + 1), prefx, statsout, region1['nsres'], masterlist);
+ #Since we are now done with the loop, close the stats file.
+ f.close()
+ grass.message('\nIterations complete!\n\nDone with everything')
+
+
+
+
Property changes on: grass-addons/grass6/raster/r.landscape.evol/r.landscape.evol
___________________________________________________________________
Added: svn:executable
+ *
Added: grass-addons/grass6/raster/r.landscape.evol/r.landscape.evol.html
===================================================================
--- grass-addons/grass6/raster/r.landscape.evol/r.landscape.evol.html (rev 0)
+++ grass-addons/grass6/raster/r.landscape.evol/r.landscape.evol.html 2014-10-01 17:24:19 UTC (rev 62151)
@@ -0,0 +1,371 @@
+<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.0 Transitional//EN">
+<html>
+<head>
+ <meta http-equiv="content-type" content="text/html; charset=utf-8">
+ <title></title>
+ <meta name="generator" content="LibreOffice 4.2.6.3 (Linux)">
+ <meta name="created" content="0;0">
+ <meta name="changedby" content="Isaac Ullah">
+ <meta name="changed" content="20141001;101831503148855">
+ <meta name="CHANGEDBY" content="Isaac Ullah">
+ <style type="text/css">
+ <!--
+ h2.cjk { font-family: "Droid Sans Fallback" }
+ h2.ctl { font-family: "FreeSans" }
+ tt.cjk { font-family: "Droid Sans Fallback", monospace }
+ -->
+ </style>
+</head>
+<body lang="en-US" dir="ltr" style="background: transparent">
+<h2 class="western">DESCRIPTION</h2>
+<p><em>r.landscape.evol</em> takes as input a raster digital
+elevation model (DEM) of surface topography and an input raster DEM
+of bedrock elevations, as well as several environmental variables,
+and computes the net change in elevation due to erosion and
+deposition on the hill-slopes using the USPED equation, and in the
+stream channels using a process equation based on the excess stream
+power or shear stress. The module has the ability to run recursively,
+looping over several iterations. The time interval represented by
+each iteration is determined by the scale of the input environmental
+variables, and as such, all input variables should be on the same
+time scale. The script creates a new map where each raster cell
+carries a numerical value, which represents the simulated meters of
+erosion or deposition (ED) estimated for that cell, under the
+specified conditions of rainfall intensity, soil erodability, water
+flow, and vegetation cover. This map of net ED is then added to (for
+deposition) or subtracted from (for erosion) the topography map of
+the previous time step, to create a new topography map (i.e., as a
+DEM) after a cycle of landuse and landscape change.</p>
+<p><b>R</b>, <b>K</b>, and <b>C</b> are environmental factors in the
+USPED equation that relate to the intensity of yearly rainfall, the
+erodability of soil, and the degree to which vegetation cover
+prevents erosion (See below for a detailed description of these
+factors). These factors largely determine the amount of erosion or
+deposition that occur on the hill-slopes. <b>cutoff1</b>, <b>cutoff2,
+</b><span style="font-weight: normal">and </span><b>cutoff3</b> are
+values of flow accumulation (amount of upslope area in square meters)
+that determine where surface processes change from soil-creep to
+laminar overland flow (sheetwash), from laminar overland flow to
+channelized overland flow (rills/gullies), and from channelized
+overland flow to full stream flow respectively. <b>kappa</b> is the
+rate of diffusion for soil-creep in meters per 1000 years. <b>sdensity</b>
+is the density of the soil in grams per cubic centimeters. <b>rain</b>
+<span style="font-weight: normal">is the total annual precipitation
+measured in meters (or the average annual rainfall in meters per
+year). </span><b>raindays</b> <span style="font-weight: normal">is
+the total number of days on which it rained in one year (or an
+average value of days per year). </span><b>infilt</b> <span style="font-weight: normal">is
+the proportion of rainfall that infiltrates into the soil and thus
+does not contribute to runoff (values are between 0 and 1). </span><b>Kt</b>
+<span style="font-weight: normal">is the stream transport efficiency
+variable that describes the cohesivness of the stream channel beds
+(0.001 for normal gravel/sandy/silt channel bed to 0.000001 for a
+bedrock channel bed). </span><b>loadexp</b> <span style="font-weight: normal">is
+the stream transport type variable that determines the type of stream
+transport modeled (1.5 for bedload transport, or 2.5 for suspended
+load transport). </span><b>alpha</b> <span style="font-weight: normal">is
+the critical slope threshold above which the model will simulate the
+cumulative effects of mass wasting (landsliding). These</span>
+measures all need to be determined empirically for a given landscape
+under a given climatic condition, but the defaults are average values
+for the Circum-Mediterranean Basin.
+</p>
+<p>By default, <em>r.watershed</em> is used to calculate flow
+accumulation modeling using the MFD alglrithm included in GRASS 6.4
+and higher. This can be made backwards compatable by checking the -f
+flag, which will use <i>r.terraflow </i><span style="font-variant: normal"><span style="font-style: normal">to
+compute a flow accumulation model using the SFD algorithm. This will,
+however, produce much less accurate results, and users are therefore
+encouraged to used GRASS 6.4 or higher.</span></span></p>
+<p>The user may use the <b>statsout</b> option to define the name of
+the file that contains the statistics of erosion, deposition, and
+soil depths over all iterations. The default name is
+<tt class="western">"mapset"_"prefix"_lsevol_stats.txt</tt>
+(in the users home directory).
+</p>
+<h2 class="western">CALCULATING SURFACE EROSION AND DEPOSITION</h2>
+<p>Because physical laws that govern the flow of water across
+landscapes and its ability to erode, entrain, transport, and deposit
+sediments can be expressed in mathematical form, they can be
+translated into a scripting algorithm that modifies raster landscapes
+(i.e., in the GIS) in ways analogous to the ways in which real-world
+landscapes change. There are various mathematical expressions of the
+relevant surface processes in the geomorphological literature
+depending for example on the processes selected to be represented,
+the simplicity of representation desired, and the degree of
+resolution desired (Clevis, et al. 2006; Degani, et al. 1979; Mitas
+and Mitasova 1998; Mitasova, Hofierka, et al. 1996; Mitasova and
+Mitas 2001a, b; Peeters, et al. 2006; Singh and Phadke 2006; Warren,
+et al. 2005; Wischmeier, et al. 1971; Wischmeier and Smith 1978). We
+use the Unit Stream Power Erosion-Deposition (USPED) equation,
+derived in part from the widely-used Revised Universal Soil Loss
+Equation (RUSLE) (American Society of Agricultural Engineers 2003;
+Degani, et al. 1979; Mitasova, et al. 2001; Mitasova, Mitas, et al.
+1996; Mitasova, et al. 2004; Singh and Phadke 2006; Warren, et al.
+2005; Wischmeier 1976; Wischmeier, et al. 1971; Wischmeier and Smith
+1978), to calculate net erosion and deposiiton across each landscape
+cell above the flow accumualtion breakpoint <b>cutoff3</b>. USPED was
+developed for hillslopes, small watersheds, and small channels (i.e.,
+rills and gullies) (Warren, et al. 2005), and is less applicable to
+larger streams and rivers. Therefore we use a different process
+equation to model erosion and deposition in stream channels (see
+below).
+</p>
+<p>Net erosion and deposition rates on hillslopes are computed from
+the change in sediment flow across cells of a DEM that have flow
+accumulation values less than <b>cutoff3</b>. We approximate sediment
+flow rate from sediment transport capacity, assuming that water
+flowing over landscapes normally carries sediment at capacity.
+Transport capacity is calculated by combining a rainfall coefficient
+(R, MJ mm/ha h yr), soil erodability coefficient (K, Mg ha h/ha MJ
+mm), and coefficient for the ability of vegetation to prevent erosion
+(C, unitless) from RUSLE with with an estimate of topographically
+driven stream power as shown in equation (1)</p>
+<p><img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAJ4AAAAWCAMAAAAGuNX8AAAC61BMVEUzMjMzMkgzM1IyNmc0O3Q1SVsySmwyR3RGMjRWNTRjNzVxOTJMMkRLNlBVM0RWOVRkOEVgQzJ4QzN+VDVMSlRSRkNQS1VKSGxBQ3VTSnJSUmRQXHJGZH1eZmxkTkNrV0JnXGN8WndsZFJ0a15obm1ocHB8cmFzenoyR4kyTZY2WIkyWas8YJs4a681Zro3bMQ7d8pZWJlcXKFLbI1beYVbe5dDdahYaa5UabpsWo95Z5hhZaZ9faVDfMdUgah/ho5qh6tsgbt5gad0nbtFgctHiNNNldtYjNVUnN5XqOVmhcVkidFsmslkldl7os18qdphs+prs+lqu+5zvO10v/F4yPSFRjKIUjKSSzKZUzOWYzipWjO0XjKraz62ZjaCWliaWF2HVGWQVGCMbFCHdVeccUSRc1yacVONZG6HYXCDYH2Ifm6jdUO4Z1W3e0O/d1+kZXCrcGm1bmu4cGzDbDfEdj2JbpOcfYKZc5qMdKelb4aBhVSEgXixjV62imC3mXG6oHfHkE7ShkXWjVvalEzZmlXKjm/Ti2LWmmbQkXXao2bgnlHlqFjipGTqsmLrs2zuumvxvWfzwGb0yHaHiIeVioiRkY2DjaqJkbCtkoK3k4WxlZu6mamxn7yJnc2Fq8iBr9GGs9+UscmUt9mJveewp8a3rtirw6yrzbyBxtSNzNeCy/SI0/eaz/SZ2fie4Pqkxd6nxOCg0O+n1vO1zuin4/q35PWw5fu06/y85Pm87v3Ago3UqoTZuIXApaPLq7PhrITgv5XdxYDswITpwpfzy4L0z5n41Yv42pfmyKngybb016fy2bv74J/746X64rb447v86rTaxN3c0uTG5/nF8/7T7vvb9uzX9/bT9P7V+P/a9f7a+f/o7NXm9NL348L55c
L55Mv97sT87Mj46tP988n+9dP+9Nv/+db++dvo+eLg8f3l+Pbk+//s/f/99eP/+uP+/ez38PH0/v7+/fP+/v78A/sAAAD////BV6T/AAAA93RSTlP///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////8A7mY2kwAABAJJREFUeJztlmtUVFUUx+/MRZkhU0HTyrTIjEdJSiYKPjIBtTQdXeCMDzRfuXpQsjTRzAeKD1CSJYK9dFmW75i1goCyp8KAyOgoUCOUQAwwGdMd79nns+c+5s65jM431/jBPV9m79nnnN/+733PHQbf18YEGsC/PcDzb//f9OfeUzy4svr4klObpwEV2v/oD3RK5xfdqiXc53YMLatPvGmAe4sHvxZjR7j5vSO/z+KpaO2Eph5JarPOBux4y5yTQuFxpl7PhGpHPK0543NIbSoz5tnHj3oEgEsLmEn1KJt5qUTRBFo/mVz4tl6rPsnarxhb4v5+o3T36bvX4NjQLe36DjNm01Ovkzq4jXZcHftHRhn24tkO805TIo9+tPuKkBXVDLnaMo+PsiLt+J/133kbhnIGC55t/HV6XYdJUw775tXG/fv+94DvZpWzpN8gO/I6rtTNBQzbU2CfwRHfSKnX2oStfVMwusr7bOA0TQfYq/1W8Y2J/IVDVBloOyvKxn1DL0Y7PwyyuzPKPzY4wqVS0Ln
CFdrS9i2PNUDN8ok/hbJnRRYpnVuYxGOLbg7JzJ9GllUl1FPqEfy9QQ3ilwOMbKzkkx5BbViiy3Outc/JT4/SIHW66b7qwPmZ2VHqVlgTbqI9xZ3b2AZoN/b+6i/TSBd2vytPk7XvGWhZ8LywwjLO5V3kweMWxkjRywcLC0Qrapaw2cJhL5oJT+dnFYKfq3nE4MERQmidplzZTcohZov/M+u5blW8+uFjfGsjziNluxeN4pGA71wqDQ3kaSOe1MwWCepim33xbCFz7zgh3OLkLjKVQoH6GSCUMTI9uN4jCAl1GSmVpByStaYYrY3xqlBN4s7FTDRpc76MB7mRzdi5TMJzpyXxXDorOtbR1KMt40EeK2sgqCd+isTibf1Pox2RAgDcEk626Q2OsInywULIEUapJOVg2D2kYFPoKF4dd58L1aSAiJcm4nnVs4WQGbTokoXV1ljqJpTx3Gs9h1QqsydoBPuDyHPeR7kaSBfKoVI3z6u0Oy1IFBNd8140VVNcuCM1Wd2OlnrcZnyBV6knzx7kC3f1AUZ4MsjsUZMt49lC7jDfInZMN+4yKkJwi0gZaB1b5k2p048l3eB+rlAibR8QtR0PzVDvWPWqC7KSgMIj74dscaTIKS64kio1KZ8mEfHQxZXM2BLfSwW3b9WxZh52aKbaPX60mSc3c/CXSgfg0qqBQyOKGhW/ZkB0Bdz4iIkuUfFZdE8UvGK/sZKZ+d9v+mBzzXzNERe2vAbirr0jhg96WTwDdtGvhgD/JXBm9nwRdGXSL70A48GFs+oRgKpjdCDQf6jcvzSp/I6vXbQbaDyMeH9uwPH8232Odxs55M9JfrWQWgAAAABJRU5ErkJggg==" name="Object2" align="absmiddle" hspace="8" width="157" height="21" border="0"></p>
+<p>where <i>A</i> is the upslope contributing area (a measure of
+water flowing through a cell) and <em>B</em> is the slope of the
+cell. The exponents <em>m</em> and <em>n</em> are empirically derived
+and vary for water flowing over nearly level ground, on hillslopes,
+in water catchments at the heads of gullies, or in small channels.
+The sediment flow rate is largely determined by the amount of water
+flowing (contributing area), its velocity (a function of slope), the
+erodability of the substrate (K factor), and the ability of the
+vegetation cover to prevent erosion (C factor).</p>
+<p>Implementing the USPED algorithm in a GRASS script combines GIS
+modules for calculating slope, aspect, and flow accumulation (the
+amount of water that flows across each cell) using map algebra. Data
+used by the script includes a map of initial surface topography (a
+raster DEM), soil erodability (a constant for uniform soil or a
+raster map for variable soil), vegetation cover (a constant or raster
+map), and rainfall intensity (a constant only). We also create an
+underlying bedrock topography map (a raster DEM) to limit the total
+depth of unconsolidated sediment that can be eroded. Soil
+erodability, vegetation cover, and rainfall are expressed as the
+K-factor <i>(K),</i> C-factor (<i>C</i><span style="font-variant: normal"><span style="font-style: normal">)</span></span>,
+and R-factor (<i>R</i>)<span style="font-variant: normal"> </span><span style="font-variant: normal"><span style="font-style: normal">components</span></span>
+of the RUSLE and have been calculated empirically for a variety of
+setting (Boellstorff and Benito 2005; MartÃnez-Casasnovas, 2000;
+Essa 2004; Hammad, et al. 2004; Renard, et al. 1997; Renard and
+Freimund 1994).
+</p>
+<p>For areas of the DEM that have flow accumulation values greater
+than <b>cutoff3 </b><span style="font-weight: normal">(ie. areas that
+are proper streams), we use a case of the transport limited process
+law that is formulated for water flowing in stream channels (Howard
+1980; Tucker and Hancock 2010). This is done by first calculating the
+reach average shear stress (</span><font face="Times New Roman, serif"><span style="font-weight: normal">τ</span></font><span style="font-weight: normal">),
+here estimated for a cellular landscape simply as:</span></p>
+<p><img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAHgAAAAXCAMAAAD3GHtvAAAC01BMVEUyMjIyMjw7MjI/Mj40MkUyMlUyMmIyNGk0O2IyOWoyO3QzPnoyQVkySWkzQ31HMjVTMjJVMjpbMzJmNTV2PDJHMkFMNVpaNUdZOlxkN0d3QTJDQ0NDQ0hPSkZQS1VHU2RNVX1nVUpkUWRnXGN0ZF50a15scHp3bmRycHp9eHB4fXozSYgySpA7V4QyUZQ0XJw/UZoyVKEzW6o6YJ0+aK81Z7c8d8hEV5xFXKlebYVXbZFfd4BbeYVEcL5ZYaNga4Jsap9ueYN0faF+d6B4ebh5gX9Mhb5ZgL59iI51hpV7lJ9qh6t4lLtHidNMktlYicxVlthSnN5WpeRdsutulcZ4nMNkot57os10rNx6r9hlo+BmquNpruVisultve9pv/FxtehtxPR6yvSERDKDRDuCTDKJRjKJTDKFVTqZTzKXUzqSXzeRWzqQZDqlVzKxXDKrYT+1ZzmQW1yDaEqGZVGAenOZdWqsY0WjdUO0cEy7c0y6ekOifnTIeD2PcoCRf42GgXeAgHuPiX2SiXm/gVmrhHW0lWq6oHfNlVXQh0fYkk3ZmFTPjWjEk2vNnnjWqX/kplrmqljorlvrsl3lrGjstmXps23uuGLvvGrjtnrxvGTzxG30yXeIiIiKiJ6Wi4aUlo2bm5Wfn5+Hj66Glbifn6OfqrOxlZuzq5WpsbO0tLWPrs6FrNiAz/aH0/eV1vaU2vmb3vq7w9Cjyemt2vW12POg4fuq5Puu6P207P275Pm86/y+8P7DoYXDpprAr6TqwobzzIX304j11p/52JP53Jnhy7zy16vx2L/746T85qv647T447v96rX87Lz/8bvB2vLY7NnC6PrD8v7M8/7Q7/3V9f7U+P/Y9v7b+P/35Mv55sP868L87sj87NH+8sX+8s
7+9dH+9dr/+NX/+d3i9v7j+f7r/P/++uT+/Or28fL0/f7//vX+/v78A/sAAAD////VKWCYAAAA73RSTlP/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////AOc7XCIAAAMlSURBVHic7ZX5N1RhGMcvpjlFoaSkjdKGtK9H+76o0HrIQYtSKZXSnna0KXTUSIsKlVYkp20iIVQzTZw4TGHu7d7n/Rt637vNjV916pw8v9zvfeZ9n8/7PssdCv0lo9rAbeD/DWw8PsttqQk/9889Ml2PEHduYNQKlvxCv7z9uFFcJerSLhQ1sqlVwMyeQPgevgSYvQEA5yeZIHlcNb07ABCqOfgUpFWi5k7Pv3ErR/LC+2gL78P9hVMKjk0W3jGDpjwBM5h+4da7J7FeywEprbDbMwQXOhQY1wYBujK82rjeDyDJRoeM2x+YsyLq8pkm5V6I9zQgrV2QHJE74FmF4LVVMMjg0n7zMzO6eGZmZX38HXxZnYeQ1iqUibR9Rm6NJX53vM/FDT3qOvGzEE/UEEd1XKZINB3hAxjsL0f8EebDEry1TgYDdhTZLxKWaCwpwdQ6pOFBdiGozLl9zD0WLqkeIlRi71/XdXZjTZjXV7Ke6FqsuXcZO6hRBhlc2j0bPmxWOAqt7hNWgipbBpP3JItQQZWlpQp2swrVOQ8zMKfwSnhjr1rJQkK7XAL20zqG4h0
qEgiJmuzmUux8AX2/85wEvKhyd1WtIiUWHXE2eh4sgCQwF6vORy0Minf2nbfFWgdvp36KpgJxsEckN/5XxJuTHIl6Nd9fsbiri3jNRPqw9L72ueSqvKMhYiwrgPKUYGarh6HZjdOFgtV09YGG8BCgI211WpItrWM2X4JKhyC+NI5mjTTjTQh+EkC5A18lnAHsIIcq6cz3Wf2a0U1KcMPGEWJjKGvM946XHlX2wHUpccqvd/YHuGCjr3DwZVGpE390UefyGUqWxgIukR7SWMq9BQnqAvKI55NgBtMRHuY2UCY7xQUXiN6GT66dXAWJ45vwHLNcvIcBEhezzLE5JklXHMqB4hlSECYST33ZOk+9FIjeOgZPG33GJUcxx3zpqWBAzc34Ki29ir/WrnnXr+LxYc5OG0K+CcyZCdfSG3F4nDdR12yg+i6UuLUn7Tq5D+iz4IsU/tsJy06pqVGD70qzLn+5uJ9sS/AftH/1T6IN3Hr2C5PS9jGdWNOOAAAAAElFTkSuQmCC" name="Object1" align="absmiddle" hspace="8" width="119" height="22" border="0"></p>
+<p><span style="font-weight: normal">Where: </span><i><span style="font-weight: normal">9806.65</span></i>
+<span style="font-weight: normal">is a constant related to the
+gravitational acceleration of water, </span><i><span style="font-weight: normal">B</span></i>
+<span style="font-weight: normal">is the slope of the cell in
+degrees, and </span><i><span style="font-weight: normal">D</span></i>
+<span style="font-weight: normal">is the instantaneous depth of
+flowing water in the cell. </span><i><span style="font-weight: normal">D
+</span></i><span style="font-variant: normal"><span style="font-style: normal"><span style="font-weight: normal">is</span></span></span>
+<span style="font-weight: normal">here assumed to be roughly
+equivalent to the depth of flow during the average minute of
+rainfall, calculated by:</span></p>
+<p><img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAIoAAAArCAMAAACQGeAdAAAC61BMVEUyMjIyMjs8MjIyMkczMlEyMlgyMmAyNGk0OWM0P2k7M2A1PXkyRWoyQXVFMjRSMjNgPDJxOTJJMkNZO1REPnFQOmBkOEVgPEh2P0FyQjN+VDVXV1dNVX1SUmRQXHJdd39kTkNoWVd6XE58WndscHp3bmR8cmF3dHUyRIYyTZQ1WY4yWqwyXLE0ZLdEV5xcV5lDW6RcXKFLbI1YY4xfd4BZZ6tUablcfqJga4JiY5JueYN6a5pgabB9a6VAfctsgbt4g790nbtHgslGidRNlNtYjNVSnchRnN5dkdZdmt1YqOVdsutkidFslc9kldlxnMV3odp9o9B+o91jreZptet3tulzw/J7y/WGRjKCUDKOVDSRSzOTVDKSXzebUDKsWzW6ajeCWliaWF2HVGWVV2KHYXCAenObYnGsY0W4Z1W/d1+kZXCrcGm1b2+1b3Cye2TDbDfLez6JbpOZc5qMdKeJeaelb4aBhVSCgXuQn22xjV6ujGa1i2O+hny0lWq4nXPLhljRhEPUiEXWjVvak0zVnlLcllzdnFDbnlrMg2LKjm/Im2vLlHbQh2TWj2HamGXOq2jWqX/lqVjoq1nrsl3ipmXst2nxvm3tyXf1yneIiYeVioiHj66KkrCxlZu0naa6mamxn7ymsbmKrMmGs9+Nv+ywp8a8w5qrw6yNzNebyNeBy/WG0veS2fma1fib3Pmf4PujyN6k08ql3Pmr1PC0yOij4vuo4vSq4vqt6P235PW27f285Pm67/6+8P7Lq7PNvazRu6DhrITpuYTdxYDqwYftwpLwx431y4H0z5n304n105L41pb52pT53Zvpyajv1Ljx1K/02av53KP74J/74qP75Kn65rfW1tbD6PrI7PzD8v7J9P/d7PbX9/bU9P7b+f
/u4Njm9NL66cb77NH/9Mb+9Mv+9NL99Nv/+NX++dvo+eLg8f3j+v/r/f799OT/++T+/Ov38PH0/f7+/fT+/v78A/sAAAD////yQWLFAAAA93RSTlP///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////8A7mY2kwAABEZJREFUeJztmHlcFGUYx0FXJU0wTJO8UtNCJI/ELBW1vELDAyW8IO/USslMvPLCqLxR16Sy1LSSTSMRDyQlj3AVPFf4eMQatYsj4I7MM3877/vOO8zS7OzI7K790fPH8+zxm3e/87zPPu8+68f/Z8xPg+Z+mSdERKAi04DC7Cz3hMjxXbHka4gCvx90T6JFdHUUSF5+5X7NKLaP8P3C5bH+rTu0/B5qIrr+TLqQkI9v8dRXWUHjLM0oedFkYVgZaoUUQ1ZNRKXf3uB5bv3bPPWS2SfUzdaKAqtiyQMmfjDAZsPPNRbx2waC5EXjvljQgGbJLQo7U1zX0vgAXAl5vfLRRWxu2yHo4809KiQv3sKZt1aH0jp2i8JMziJXpT6R9ly3DOVScS0q3XFcqKRmeFOudbUSj19EVtjrzxVh9KvnHmXqERwdCQOYCU63e6FNc2IvFzNTsl2ICoJigDc/fRM9tnQuIz4/aDSmdbx3kJsbTuWas1LY1MQlv3jrUUXwAHhuZTjeEwvOiqVrOTzAHw8bWiU
lNYrUjMLNNOHLtgl3lh9sUhaxqiJmPEmCuUel5JH6fN8K/u74AXQ73aN8htdh3wkvF9aMVKxaHj5VE11r8ushQKBvSh5ZyYdCwf7zVJRmFN48TND+vSyw4T6A5Nr9iwun/TjppyUDQVUEclVe8HBUmtw6k+QFpIvPdjwOJYv8OmaCVhRmntXpua1dxqzdp6MrVUVKKqL5l1IyDWfQ2f1OGTBH/DXl8Pp0dZGCCs7tAerVUODSRP/uSe17ZirUAntCfpjC1jGWiHuzj4GaSFDFFvR2Vtn3VkheDYXnEkOtQvMLHqPAzMn52PePGkfb2mVW1zmJBNUxY4yzigg45bqXo9yNe0MQsXMb3HQp9bZRlPxgtK+QWluo77xafsTq+ZSL1koyPiBhs79wWhR+lUbsa1fF7k0UJh73JW5FnaPV3veJyVEsQbGowuxx6LSQsvJNme9RILXebyikBKBjePvjrBVHAsqGY02LI676j49QoOSTWk+mpU1/YbfL9uMrFE+afbH4exHON8mpitym56ePdN3ddKIoDhrsybGdCIp9XJ1sKYLxJavj3Rj13deRFcVBg0sMwyjc2oX1s6Voj4sCMLrp5DpQFAcNLhFnBf4Y8SVCEaM5UOic5sB05YX0oygOGiJKyXzrdoRCIhgRsCVoqOoO6aiVaoMGmSgICrv0MI9QxAiphl8QSpSXUOSDRumcdDJm8IBQ4JRQzAiBxi2GA2gO8VZW5IMG+0EOHjPErDgSSLceJMbB5kCTF2vFadCw9aEjHkaBO0VFRRsN+27TaAsZKlS4qylKL0rVoAEXdy0eQlNPv8w8rhUawdilzGt9RTZoXO9XsZX+UcHmhgSIx4cchWfXvubFbkuN+9zEzMjRvYwnUJjJOWcibnvgRPdAVpa/8sOrN/STePxk1mH/oyjZQ0WL90bEF7SYAAAAAElFTkSuQmCC" name="Object3" align="absmiddle" hspac
e="8" width="137" height="42" border="0"></p>
+<p><span style="font-weight: normal">Where: </span><i><span style="font-weight: normal">R</span></i><sub><i><span style="font-weight: normal">m</span></i></sub>
+<span style="font-weight: normal">is the total annual precipitation
+in meters, </span><i><span style="font-weight: normal">i</span></i>
+<span style="font-weight: normal">is the proportion of rainfall that
+infiltrates rather than runs off, </span><i><span style="font-weight: normal">A</span></i><span style="font-variant: normal">
+</span><span style="font-variant: normal"><span style="font-style: normal"><span style="font-weight: normal">is
+the uplsope accumulated area per unit contour width at the cell, </span></span></span><i><span style="font-weight: normal">R</span></i><sub><i><span style="font-weight: normal">d</span></i></sub><span style="font-variant: normal">
+</span><span style="font-variant: normal"><span style="font-style: normal"><span style="font-weight: normal">is
+the number of days on which it rained in a one year period, and </span></span></span><i><span style="font-weight: normal">1440</span></i><span style="font-variant: normal">
+</span><span style="font-variant: normal"><span style="font-style: normal"><span style="font-weight: normal">is
+a constant relating to the number of minutes in a day.</span></span></span></p>
+<p style="font-variant: normal; font-style: normal; font-weight: normal">
+Then the transport capacity is calculated by:</p>
+<p style="font-variant: normal"><img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAE0AAAAdCAMAAADYWRF1AAACZ1BMVEUyMjIyMj07MjIyMkU6MkMyNlEzM2IyM2gyO2wyO3Q1SVsyR29DMjJFMj9SMjJUMj5cMzRoNjJQMkJgQzJ/RDJDQ0NDQ0hPSkZSRkNkUWRsZFJwYldobm1ocHAyQ4IyRokyTIMyUIYyVKEyXK44a684bb9VW5tbe5dDdaJ5Y5djaKJBesFUgah/kqJCgc9DhdFMktlUouJkgsxtnMR/l9Zkot5rtOpouu1wtehzueF1vOx0v/F6y/WUTDKcUTKTYDeVYjiiVDKrXDKxXDKwYDO6aDaKU1CccUSacVONZG6DYH20fEO2bWi4cGzCdz6BZJKcfYK6n3LQgkLUjV7alE3XnVnao2bkp1jkp2XlqmnqsmLuuWT1yneHh4eKioqWjIqHj66JkbC3k4WvirSsi7q0rpGpsbO0tLWIi8yJnc2BktOFnNiWiMSArs6Ept6UscmTuN+GqeCEu+y+tdyrzbyBxtSBz/aPyvWVzvWlwt6vyOCm0fGk2vmv1fWr3/qxz+27zeSz0fK73Pet5/217P275fu/6f3Ago3MlaXIlKnDoYXXqa/draHSp7nhtITnt5bhsp7Cq9HFv+HlwpzxyIHxxpL0ypP205n52JL63prgybbyzKP11Kv63KTx2bz74KH956v24b/75LT96rT967n/8L/SwNzaxN3B2vLF4/fC7v7O6/nD8v7N9P/V7/nb9uzU+P/a9P3b+f/u2Nfj3u/55sP44sr86sP868n56dD669z+8cP+9Mv+9NH/+NT/+d3l+v/p9f3q/P/46eL/++T//er29fX1/v7+/vT+/v78A/sAAAD///+WEM3BAAAAy3RSTlP//////////////////////////////////
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////AL7mpKIAAAJDSURBVHicYzhFTcAwahoucHLzIRT+to2nTm71LPEwO0mOaUcmnUTh7208eWqPV2WYNVmmba9F5R+dePzUDLX1Ps2nyDFtdhsq/0T3ppMRZrs1N5DjtpM9G9AE+lcd9W2er7WSHLed6NyEJjKtBsEm1bT9KYfQRKaXk2/avpzjaCJzixCRDDHtiB2XjDCrrBRLGbruk/McGHRWnghlVK8Ha9qTgW7aojx00xZnHt9vp3/8RPtaDMecCJLbcGpHQDVUy9IsoGlHuqTFRUFAzPTkqSXZx9FMAybo+fw2p04sR7cYGFC2+scnJ8IjcrEJ0NgFkrppqQIKaenpK06eWmqCbhrQR+Gca8CMGCYGCGCH8GcKlSYVINQvBuk9CcQz+Qwgzt2dgWnaESdVSGQtio+LBYOETRBbmCXMkPLSEqhLTkYz20BEdmdhmrZYwBw1/8FscVTxlluF4MNccjKIrQVqfjZ6LJw6GQWTBLoNAhIPgLgLBMy3i2gj0hgsTo84y0HDcnEhhmlHXZUPQFio4Qay5eQUXiuwhhNbTsLT2347JagGzPR2aomgIXaPugBtOerHASohTvYB9R1JPgANAXmo26ZXIJSDTTsxy51BowEzdZzaGcijWHV8oQMDd/6BU9vdcg+dOtEByacnQxgtIPZPbUYzjTiwxx/ompMToAn8xLFjYNPgAiSaNsMI5Po5raiiJ3oPkGPayUjLgyAXNqIKHylGCiHiTTsRbAxyxY
llqOG7owmJQ3qddRjVtF3I5d1grp0BWs1h0gXeQP0AAAAASUVORK5CYII=" name="Object4" align="absmiddle" hspace="8" width="76" height="28" border="0"></p>
+<p><span style="font-weight: normal">Where: </span><i><span style="font-weight: normal">K</span></i><sub><i><span style="font-weight: normal">t</span></i></sub><span style="font-variant: normal">
+</span><span style="font-variant: normal"><span style="font-style: normal"><span style="font-weight: normal">is
+the transport efficiency factor related to the character of the
+stream bed (0.001 for normal sediment to 0.000001 for bedrock), and </span></span></span><i><span style="font-weight: normal">n</span></i><span style="font-variant: normal">
+</span><span style="font-variant: normal"><span style="font-style: normal"><span style="font-weight: normal">is
+an empirically determined exponent related to the dominant type of
+transport in the stream system (1.5 for bedload transport or 2.5
+suspended load transport).</span></span></span></p>
+<p>Net erosion and deposition rates are then computed across the
+entire DEM as change in sediment flow in the x and y directions
+across a cell as follows”</p>
+<p><img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAM0AAAAnCAMAAACbvTy/AAAC5VBMVEUyMjI2MkY0M1Y0NWYyPHc0RFoyRWoySXY4U3FFMjRWNTZjNDNJM0hJNlJYM0RZNlBINmBjOU9WQDJiRTN9SzhcQ15XVEhVVldLV2pNVX1XTHRGZH1rWFR7QkNsZFJ6aFlocHB3bmR2Y314d3cyR4kyS5E/X4QyVJcyWqs4Y6g3abo/dLk3bMQ7dslHRoFAUY9EWp1cV5lcXKFIdoxXYIhfd4BEb6lDdaJGc7VTZqxUabpfeaVga4JqYqZ1bqN9faVNfMR3holxhJl5g6l1ma1Egs5JitRLkNdYjtdcmtpWpuRfsupmiMRkidFsn8Rhm9R+o91hp+NouexrvvB3tudzwON3xvOJRzaMUjKSTjOZUzKdZz2pWTO0XjKpZTq5ZzaTWEKaWF2QVGCLYUSMbFCHdVefa0SccUSXcFODeGWId3inZlOvfli1bkq4Z1XHeD7DdEGVaoKZc5qMdKelb4agf62BgHu4g1atiX27hGC0lWqyknG7oXjHkE7ShUfWjVvak0zcmlrIm2vTi2LQkXXfr2PgnlHlqVjmrWXtuGjpu3/xvWXtwX3zwGb0yXeHiIeHk5uVioifnpqHj66IjLOKkrCVgbefn6aVobqriYG7hZ21l5mznKW/mqe6mamxn7yzrJW7oqWJnc2Cn9eBr9GVptiXuNiJveeUueeolcijo8WypsiBy/SI0/eew+uaz/SW2fee4Pqz1Mmjx+Sp2PS2y+O02vSo5Pu46fu88P7Ago3Wn43XqIXZuIXds5vAo6XLq7PKs6XUoqLhq4XhupXgvqLqwobtwpL1y4H0z5n31Ir31pTix6bgybbo0rX32qTz2LP74J/64qb75rb+8r7NwtPaxN3b2OjL5ujD5frG8/7W7PjY9/7o48Ho7NX658j26t
T+883+9dT+9dv/+NX/+dvi4eno5/Tg8f3j+/7v9/Hu8vnr+/799uT++uP+/Ov0/v789fH+/fT+/v78A/sAAAD///9mta+aAAAA9XRSTlP/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////AMLtd88AAAV8SURBVHic3Zn7W1RFGMd3WUgXMQwqyqyopAuBlRWlXc0sFbJ4MDaQwm6mVkZWVFaElWyEttXT/W4lWEJiorRd7ILVthjwLFvGurWysMXmcs7hzPwBzczeZnfPXmYZftD3B2DgzGfOd+Z933n3RQWPJFNNMH90JNqAF5QeMaoZOwSYnh/6wE2NhPejyYnJBY2zjNGhNJVJjfCoVlXEsr+ga0fIuPcmKQku6L3KHgNKUVnUgK5mt+O+UjRVqEzLzU7JPUXdFnvG4EMhuwiFh23j54ZDKSqLGvlf5A6uRxCr/yPJVblAEr9SejnKTItCHQg83Tp+bjiUojJnAflHdOqOEWhJL4fyH36u8y3DjBIPHHrMcOpldih3Ge5O6cTLlPvWGwNQ/Adt/qaF0eJDkesFHVw33QZ67rr62yxNO1CABqlMnta3+qVX3txMfBjoU6n9E+4tByZtmXh/KRCqzrFb5rjl55Bzi2vaybN7a7OL3E0aFMrmuR4WLgF1jD6lMcoO3aTmvyvzPQrQIJVBDeg+rgNAy5VkVaEq30c
YescITbNGoGvzQPdUdCAb1W3mzBbJcQB5zx2deN7GOW5L+s01C9A2WmZHBrsyF2MhAdlgk8YGxVsLJbk+z64ADVIZ1LgqrkNTXZ+SBNI/pcx3uhZtKdDn4VcBjSlG/JfFrmWqM9ELwmGy8LCuDMg1k6fjXbRcYE+Qi7FQ8IKImupCCTSgdSKhQSqDGnMGPuL+FuB98S99v0buC/T4wMGYHgeLdUopEL/JUrcC6Lq9k8xDTmfSFuN5ltnuBLk4KlDgEFBT6gA5G6ImEhqkMqghPip/TXZDrD6bfq0XMda6FYUOcpnMdscAdOgKJCivwZm2CWkEL6QTnzDPjbxwYnAxqKJA8nsaaJhpR3ETDg1SGdSY8Kq9H/scYjGdnJw69bUfnm4crTrLDp6/1GMu8YC6hUjNeuw2JnWb/FPzBpLkmsKza2yuuUSS664HWG9ADQiHUlQGNc7lN/6371Ucb/Jv96gu2UXvsvVO1eTtAA6uPOGM+XZo1p5kmI9d2XwDWsZ5W9o1u4Gz4lwPBM8qXIsxuF6QY5lqiXtv1qSd+25Rf+KJgFJUnxrQt0J1cW1ttsbYtxz9cNrlOxXuhdEf3tilWJlENdfa0KAfXhtIafL+wAqM3HAoRfWfDajPQ3ne+bpHrskbgeDnjKVs5aWygV+20xjwa0tgaL0iMh8kBaWpfjVCFQ4oVGM4dThfyjWpcYqWxEz8nt5H57bg3Zm8mjAoTfWrsRzdCvu3Io3dWuyEQK9C30wqn2mSliYfogeUP41DTSiUpvrjRp8ybxXO+6B+po2M1aga6n/Z4LXXOH3Momw8aqKaT41YXSCJb6MDdFUWklurjuS+kCe5GeyekZNzvDonJ+ciGzcurQYXFCTLoBIWqxnU5bvDzobTqv6FA2fDE+r9Cho1uKAQtkmkRoGgQYOPxqT1PZp83ES1CfQ0mRQUwuNLUQmLf3jmxA4eCTqWTZyag09q0+YZVmVpvlinPcpgWJ27
ReFTSBxj7X8kpoaRyqkDxdz/mBgqHzXBPgVPY6fyURPsU/A0dirHXifpU3A3JionTwv0KRxPTGt3rsjjEkJU9wP0rDzfDhrPi6OMixq6TyE/UPTun/u5FOB09wOapxmh48E49x4XNSF9ik2a8JqIBxVaj2mH1i1xMgIXNXSfAnZn7IjzeDJUOHhyq/hevBDioiakT/HdscV86giaik+qeE9HvClc1AT7FEMH9uyuy+8xxp/DQsUmVl34edxd4qIm0KcADWktsnkql0/hVPcDm1yzKP41yidDs/c/mKmOzxKoHSf6P4V8TPz9r0TEHCZqzJlLErqPDw81idr/U/vOq8tgJ5cAAAAASUVORK5CYII=" name="Object6" align="absmiddle" hspace="8" width="204" height="38" border="0"></p>
+<p><span style="font-weight: normal">where ED is net erosion or
+deposition rate for sediment and </span><em><font face="Times New Roman, serif"><span style="font-weight: normal">α</span></font></em>
+<span style="font-weight: normal">is the topographic aspect (i.e.,
+direction of slope) for a cell. Whether flowing water will erode or
+deposit sediment in a particular cell is determined by the </span><em><span style="font-variant: normal"><span style="font-style: normal"><u><span style="font-weight: normal">change</span></u></span></span></em>
+<span style="font-weight: normal">in sediment flow (transport
+capacity) from one cell to the next. If the transport capacity
+increases (for example, due to an increase in the steepness of the
+slope or amount of flowing water), more sediment will be entrained
+and erosion will occur; if the transport capacity decreases (for
+example, due to a decrease in slope or water flow) sediment will be
+deposited.</span></p>
+<p>The output of this GRASS implementation of these transport
+equations must be modified in several ways in order to make it
+appropriate for landscape evolution simulation. First, because of the
+way slope is calculated in <em>r.slope.aspect</em>, the flux <i>T</i>
+is actually calculated one cell downslope from where is really
+occurs. This causes problems when USPED is iterated over many cycles,
+and creates oscillating "spikes" in positive and negative
+flux values resulting in the calculation of alternating deep pits and
+high mounds at sensitive areas on the landscape. To overcome this,
+<em>r.landscape.evol</em> uses a nieghborhood algorithm in <em>r.mapcalc</em>
+to put the calculated value of <i>T</i> back into the cell that is
+most uplsope from where it is originally calculated.
+</p>
+<p>Additionally, control must be kept for the amount of erodible
+sediment available to moved. <em>r.landscape.evol</em> explicitly
+tracks this by taking the difference between the input bedrcok
+elevation DEM, and the current surface topography DEM, and creating a
+map of "soil" depth. This map tracks the amount of material
+assumed to be available for entrainment and transport by surface
+processes. A simple logical algorithm is used to prevent unduly large
+amounts of erosion from being calculated in areas devoid of erodible
+materials (ie. at bedrock outcrops). Where this condition occurs, <i>K</i>
+or <i>K</i><sub><i>t </i></sub>is made to be very small, resulting in
+only extremely small amounts of erosion.
+</p>
+<p>Another major issue is that the total flux <i>T </i>is in units of
+Tons/Ha, which means it must be converted in order to calculate the
+change in elevation at each cell (<i>m</i><sub><i>vert</i></sub>).
+This is done via a simple algorithm that uses the density of the soil
+and the cell resolution:</p>
+<p><img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAK8AAAAVCAMAAAAHHU5UAAAC6FBMVEUzMjM1MkkyMlMyNmU1Pnk1SVs1Qns9VWtGMjpSMjNmNjZFMkJUN1REO2pUOmBjO1VgRDN1QjZJQndKXW5TU2lZYGtdd39kTkNkT1B7Q0dwVGZ3WHJsZFJobm1rcHh5c2d9eHAyQ4Q4TZQ5U5k1WqoyXLE4Yqo3ZrlcV5lDW6RcXKFLbI1JepZBbb9DdaJdZKVUabhcfqJReLJiY5JueYNmeJtkZKN9a6VBechcg59Vh7Jlm7p4g797n6l5mLd9orZFgctGiNNOlNtThcRYjtdamNlWotFXpuNdsepkgsxkidFlmtZ8hsF9o9FkrOVpt+xsv/Fxr+F0ueZ5wut5yvSFRTSBVjaSTDKVUzKbUDOZUD+oWTW3ZDe2cDyKSUCMWUaKU1CfWkGaWF2XWGOccUScfFmRZG6bYXGsY0W6aFS6e0a3dVe+b2S2d2izfXvAbT/Mez7MfEONaIGBZJKOdYCZc5qMdKelb4aujGauiHKyiGa4nXPJg0fHiFjTikjVjVzZlU3cmFPMgWHNnnjTi2Len2nWmX7drGrVqH3lqVjqsl3lrWnrtmXtvXbxvm3xwW70x3WIiIePgZKVioiUn5mHj66Jn6+KkrCfqKmXpr+Sv7+rlJOxlZuhgrK6mamxn7ynoJmrq6OouqOnsrq/rLK7tLSJnc2HrsuGs9+Yud6Mq+GFv+yhosuwp8a8w5qrzLGDzPSF0feaz/SX2vif4Pumyd623cmlxOSj1/e7zeSm4Oim4/u46vvAgo3ZnYDEqYXLtZrRuZ/Lq7PNvazhq4XiuIfKt8XowY3rxpT1y4H0z5n304f53Jjuy63gybbt0LD01qb74qT86LbIw9naxN3B2uLM59jF5/nI8//V7fPU9P7a8/3b+f7h1t/j3u/66cf77d
L99Mz+9tP99Nv/+NT++dzu+eXg8f3j+v7q9f3r+/Pr/f736+b99OT/++X+/Ov38PHz/f759Pn+/fP+/v78A/sAAAD////1WLTOAAAA9nRSTlP//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////wBKPyCpAAAEYElEQVR4nO1WaXQTVRSepIhNUKSKoihWoFBQK6ioR9Ge0KqoCEVskdqIuwUUFXeRloAbKbEWcN8FN8Ql6QIqcV+gJGgzqC2daknSkdDoDHPenf++N29elpl42pOo53iO98fcO3fe3Pvd7913Zzj1vyVc7iH+2J97jEGnyAmvePVcUMVXY7kDGkCkV7p1Kye8oeF+Fb5pJWZk+dQ4VmjNiSfNVuhTZht1VplmQ1Z4pepDxoywjjnB4sM34JkcU8P3EXrlr6uK4ppHkBbN0WIz26hZqOjyowtHcTPSC4CO+dwxhSPfMZYl3d+VFV7+LUWsdijoU/I6uqNMUdtnaQjQnRMw3kjVTIxriBab2UatR0IrJsZU6dZKQwKoHycgl3W70b3SmxXeXkEN2SpV9CMBGT7OhwPRjBRvIB/fBYZpLmYbtR5JqimJq7C725BArC4DaLC2GhM3XghJvNHPj71yx4iS/jW2mWBcqLZzuuRp1IA7QVHgyC5VXuyjeJdhvOCxtOBes10CasKuYL6UZ1ocl+WKDM0cOqw
VOgrKcY/te/G5m8uV3nueml9MmNCOB+M3OHzjZlfRF23OGWa8/IZ1VJ4WKC/FcZYRMyTesF3Hi/sX3ISYkFZ1wjZqiL7k1wIt4WYlAetOcOetG3XKNuznR/vBM75nlUPh38XhQpN6knjBPXTzb85Tt4ULjO1kFt5eodck1WCu+m5keAm/2kYGKb/MbrC2UG2hGoI2eujkOu6iBGDdKS0qExc4sFe8fg50jG6GRssWbXFokpDEKzsnCry9PB443J8B4Ib1lN9nSYXgsbI1/BEYahq/uD295KK1CLMD+T6qh1GtwkEy3A5gcNfo5yq61AsHgbLhQ66xuOHaLRe89vx+bY3Wy6GTY0m8vL0SmqwtqJ5tNRP0K6hNaf0rOyfonwd4chwuAKX2rxomM6CBNjizwwXpmr7d8RG+NBHKSdClOgfQiJ8HSU2useQk4ir6FhSTTQhMVZJ4m4Z2oRXFeFJdvCftCMC3W4z9zNvZkUTLHGTbV83V51kRrgPWTo6RGRu+9owYsxPaQzWN7CpXVFQ/njZl+BydA9mJGROvcij4tHoh+sbPL8fAXUIGe+p8wJmVSFUF9C08T4CdZ8fWzvth04Mv3PLh8VuvezudcPT9TdxpH9Oa+hbSsXUZCSR/WZBH1soPTyPfMN5uaaU286VoIlINd+i0wrM6CfTdmx7Qwex7KH9Im4JcltLu6JL8o0oF3nbu6+sJ9/AI3Ubz/BVru/iNl/ffvvWutl/CdxvHo1nEe3syeOEz48T/SwmdH3/UO9CiCMtiwotu++S7uulvdganCOquSwf+3MPOZvMMVHvfG+yPAlrtFWszHPP0JLtYY5rxri7du9gLvz9egTtvXv/ACeWvzJsQ3RPPsDKjSLX+HVP2Zig5VSLvs3hmvI81wxMjn/kJH1m0cvqBQWREWf5z6W/Xnf7BmZ2DTvE3/K//q/I/3n9W/gTfyip8oJmteAAAAABJRU5ErkJggg==" name="
Object5" align="absmiddle" hspace="8" width="174" height="20" border="0"></p>
+<p>Where: <i>10000</i> is the number of meters per hectare, <i>Sd </i>is
+the density of the soil, and <i>Res </i>is the cell resolution
+(width). In order to convert the output back to Tons/Ha (standard
+rate for USPED/RUSLE equations), you can multiply the <b>netchange</b>
+output map by "(10000 x resolution x soil density)" to
+create a map of soil erosion/deposition rates across the landscape.
+</p>
+<h2 class="western">NOTE ABOUT VARIABLE CLIMATE PARAMETERS</h2>
+<p>r.landscape.evol accepts an external “climate file”, which
+should be a comma separated plain text file with four columns in the
+order of, "rain,R,storms,stormlength" (without headers).
+Each of these columns must exist, although there need not be values
+in every column (i.e., you can enter a single value for any of these
+parameters in the command line, and combine that with populated
+columns for the other values). Note that the climate file must have
+the same number of rows as there are iterations of the simulation
+(“years”).</p>
+<p><br><br>
+</p>
+<h2 class="western">SEE ALSO</h2>
+<ul>
+ <li><p style="margin-bottom: 0in">The <a href="http://medland.asu.edu/">MEDLAND</a>
+ project at Arizona State University
+ </p>
+ <li><p><a href="r.watershed.html">r.watershed</a>, <a href="r.terraflow.html">r.terraflow</a>,
+ <a href="r.mapcalc.html">r.mapcalc</a>
+ </p>
+</ul>
+<h2 class="western">REFERENCES</h2>
+<p>American Society of Agricultural Engineers 2003 Honoring the
+Universal Soil Loss Equation: Historic Landmark Dedication Pamphlet.
+Purdue University Department of Agricultural and Biological
+Engineering.
+</p>
+<p>Clevis, Q., G. E. Tucker, G. Lock, S. T. Lancaster, N. Gasparini,
+A. Desitter and R. L. Bras 2006 Geoarchaeological simulation of
+meandering river deposits and settlement distributions: A
+three-dimensional approach. Geoarchaeology 21(8):843-874.
+</p>
+<p>Degani, A., L. A. Lewis and B. B. Downing 1979 Interactive
+Computer Simulation of the Spatial Process of Soil Erosion.
+Professional Geographer 31(2):184-190.
+</p>
+<p style="margin-left: 0.5in; text-indent: -0.5in">Howard, A. D.
+1980. Thresholds in river regimes. <span style="font-variant: normal"><span style="font-style: normal">Thresholds
+in geomorphology</span></span>, 227–258.
+</p>
+<p>Mitas, L. and H. Mitasova 1998 Distributed soil erosion simulation
+for effective erosion prevention. Water Resources Research
+34(3):505-516.
+</p>
+<p>Mitasova, H., J. Hofierka, M. Zlocha and L. R. Iverson 1996
+Modelling topographic potential for erosion and deposition using GIS.
+International Journal of Geographical Information Systems
+10(5):629-641.
+</p>
+<p>Mitasova, H. and L. Mitas 2001a Modeling Physical Systems. In
+Geographic Information Systems and Environmental Modeling, edited by
+B. O. Parks, M. Crane and K. C. Clarke, pp. 189-210. Prentice Hall,
+New York. 2001b Multiscale soil erosion simulations for land use
+management. In Landscape erosion and landscape evolution modeling,
+edited by R. Harmon and W. Doe, pp. 321-347. Kluwer Academic/Plenum
+Publishers, New York.
+</p>
+<p>Mitasova, H., L. Mitas and W. M. Brown 2001 Multiscale simulation
+of land use impact on soil erosion and deposition patterns. In
+Sustaining the Global Farm. Selected Papers from the 10th
+International Soil Conservation Organization Meeting, May 1999,
+Purdue University, edited by D. E. Stott, R. H. Mohtar and G. C.
+Steinhardt, pp. 1163-1169. USDA-ARS National Soil Erosion Research
+Laboratory, Purdue.
+</p>
+<p>Mitasova, H., L. Mitas, W. M. Brown and D. Johnston 1996
+Multidimensional Soil Erosion/Deposition Modeling Part III: Process
+based erosion simulation. Geographic Modeling and Systems Laboratory,
+University of Illinois at Urban-Champaign.
+</p>
+<p>Mitasova, H., C. Thaxton, J. Hofierka, R. McLaughlin, A. Moore and
+M. L 2004 Path sampling method for modeling overland water flow,
+sediment transport and short term terrain evolution in Open Source
+GIS. In Proceedings of the XVth International Conference on
+Computational Methods in Water Resources (CMWR XV), edited by C. T.
+Miller, M. W. Farthing, V. G. Gray and G. F. Pinder, pp. 1479-1490.
+Elsevier, Chapel Hill, NC, USA.
+</p>
+<p>Peeters, I., T. Rommens, G. Verstraeten, G. Govers, A. Van
+Rompaey, J. Poesen and K. Van Oost 2006 Reconstructing ancient
+topography through erosion modelling. Geomorphology 78(3-4):250-264.
+</p>
+<p>Rawls, W. J. 1983 Estimating soil bulk denisty from particle size
+analysis and organic matter content. Soil Science 135(2):123.
+</p>
+<p>Renard, K. G., G. R. Foster, G. A. Weesies, D. K. McCool and D. C.
+Yoder 1997 Predicting soil erosion by water: a guide to conservation
+planning with the Revised Universal Soil Loss Equation (RUSLE). In
+Agriculture Handbook, pp. 1–251. vol. 703. US Department of
+Agriculture, Washington, DC.
+</p>
+<p>Renard, K. G. and J. R. Freimund 1994 Using monthly precipitation
+data to estimate the R-factor in the revised USLE. Journal of
+Hydrology 157(1-4):287-306.
+</p>
+<p>Singh, R. and V. S. Phadke 2006 Assessing soil loss by water
+erosion in Jamni River Basin, Bundelkhand region, India, adopting
+universal soil loss equation using GIS. Current Science
+90(10):1431-1435.
+</p>
+<p>Tucker, G. E. and G. R Hancock 2010 Modelling landscape evolution.
+<span style="font-variant: normal"><span style="font-style: normal">Earth
+Surface Processes and Landforms</span></span> 35(1): 28–50.
+</p>
+<p>Warren, S. D., H. Mitasova, M. G. Hohmann, S. Landsberger, F. Y.
+Iskander, T. S. Ruzycki and G. M. Senseman 2005 Validation of a 3-D
+enhancement of the Universal Soil Loss Equation for prediction of
+soil erosion and sediment deposition. Catena 64:281-296.
+</p>
+<p>Wischmeier, W. H. 1976 Use and Misuse of the Universal Soil Loss
+Equation. Journal of Soil and Water Conservation 31:5-9.
+</p>
+<p>Wischmeier, W. H., C. B. Johnson and B. V. Cross 1971 A Soil
+Erodibility Nomograph for Farmland and Construction Sites. Journal of
+Soil and Water Conservation 26:189-92.
+</p>
+<p>Wischmeier, W. H. and D. D. Smith 1978 Predicting Rainfall-Erosion
+Losses - A Guide to Conservation Planning. USDA Agriculture Handbook
+282.
+</p>
+<p><br><br>
+</p>
+<p><i>Last changed: $Date: 2009-23-1 (Fri, 21 Jan 2009) $</i></p>
+</body>
+</html>
\ No newline at end of file
More information about the grass-commit
mailing list