Simple Comparision of Two Least-Cost Path Approaches

Submitted by dylan on Tue, 2008-02-19 08:19.

Least Cost Path Comparision MapLeast Cost Path Comparision Map r.walk in blue, r.cost in orange

 
Premise:
Simple comparison between the least-cost path as computed by r.walk vs. r.cost. Default parameters, unity travel friction, and λ = 1 were used with r.walk. Slope percent was used as the accumulated cost for r.cost. Test region similar to that used in a previous navigation example.

 
Results:

  1. r.walk
    • Shorter total path (6.03 km)
    • Tendency to stay mid-slope, and cross large drainage at the gentlest possible position
  2. r.cost
    • Longer total path (7.12 km)
    • Tendency to follow drainages upstream-- not always possible in steep terrain

 
Code:

 
Prep work see attached files at bottom of page to re-create

# generate cost "friction" maps
# unity friction for r.walk
r.mapcalc "friction = 1.0"
# use slope percent as the friction for r.cost
r.slope.aspect elev=elev slope=slope format=percent

 
Compute the least-cost path via two methods

# compute cumulative cost surfaces
r.walk -k elev=elev friction=friction out=walk.cost start_points=start stop_points=end lambda=1
# Peak cost value: 14030.894105
r.cost -k in=slope out=cost start_points=start stop_points=end
# Peak cost value: 11557.934493

# compute shortest path from start to end points
r.drain in=walk.cost out=walk.drain vector_points=end
r.drain in=cost out=cost.drain vector_points=end

 
Vectorize least-cost paths for further analysis

r.to.vect -s in=walk.drain out=path_1
r.to.vect -s in=cost.drain out=path_2

# add DB column for length calculation
v.db.addcol path_1 columns="len double"
v.db.addcol path_2 columns="len double"

# how long are the two paths?
v.to.db path_1 option=length column=len units=me
# 6.03 km
v.to.db path_2 option=length column=len units=me
# 7.12 km

 
Simple map

d.rast shade
d.vect start icon=basic/circle fcol=green size=15
d.vect end icon=basic/circle fcol=red size=15
d.vect path_1 width=2 col=blue
d.vect path_2 width=2 col=orange

 
Sample the elevation raster at a regular interval

# add a vertex every 100m along each line
# this will create two tables for each
# path_1_pts_1   <-- original atts
# path_1_pts_2   <-- original cat + distance along the original line
v.to.points in=path_1 out=path_1_pts -i dmax=100
v.to.points in=path_2 out=path_2_pts -i dmax=100

# add column for elevation values
v.db.addcol map=path_1_pts layer=2 columns="elev double"
v.db.addcol map=path_2_pts layer=2 columns="elev double"

# upload elev values for profile graph
v.what.rast vector=path_1_pts raster=elev layer=2 column=elev
v.what.rast vector=path_2_pts raster=elev layer=2 column=elev

# dump profiles to CSV files
db.select path_1_pts_2 fs="," > p1.csv
db.select path_2_pts_2 fs="," > p2.csv

Least Cost Path Comparision ProfileLeast Cost Path Comparision Profile

 
Plot an elevation profile for each path

## startup R
R

## load data from previous steps
p1 <- read.csv('p1.csv')
p2 <- read.csv('p2.csv')

## can't guarantee that the two lines will be in the same direction (a r.to.vect thing?)
## p1 is correct, p2 is reversed

## make a nicer plot
## library(Cairo)
## Cairo(type='png', bg='white', width=800, height=400)

plot(elev ~ rev(along), data=p2, type='l', col='orange', lwd=2, xlab='Travel Distance (meters)', ylab='Elevation Change (meters)')
lines(elev ~ along, data=p1, col='RoyalBlue', lwd=2)
legend(0, 3000, legend=c('r.walk', 'r.cost'), col=c('RoyalBlue', 'orange'), lwd=2)

AttachmentSize
elev.ascii_.gz543.35 KB
start.txt35 bytes
end.txt34 bytes
‹ Mapping Wifi Networks with Kismet, GDAL, and GRASSupVisual Comparison of 2 Raster Images ›
( categories: )

Updates

Submitted by dylan on Thu, 2009-03-12 18:53.

Just re-compiled GRASS with Colin's recent patches to r.cost, r.drain, and r.walk. Compiles fine against develbranch_6 on linux, with visibly different results. Using the output direction update to r.cost/r.drain does not appear to differ much from the original implementation. The results from r.walk/r.drain using the output direction patch results in considerably differences-- that appear to be non-optimal in this terrain.

# re-compute with patched modules
r.walk -k elev=elev friction=friction output=walk.cost.new outdir=walk.dir start_points=start stop_points=end lambda=1
# Peak cost value: 14030.894105
r.cost -k in=slope output=cost.new outdir=cost.dir start_points=start stop_points=end
# Peak cost value: 11557.934493

# compute shortest path from start to end points
r.drain -d input=walk.cost.new indir=walk.dir output=walk.drain.new voutput=path_1_new vector_points=end
r.drain -d input=cost.new indir=cost.dir output=cost.drain.new voutput=path_2_new vector_points=end

# add DB column for length calculation
v.db.addtable path_1_new
v.db.addtable path_2_new

v.db.addcol path_1_new columns="len double"
v.db.addcol path_2_new columns="len double"

# how long are the two paths?
v.to.db path_1_new option=length column=len units=me
# 4848 meters
v.to.db path_2_new option=length column=len units=me
# 7000 meters

# map
d.rast shade
d.vect start icon=basic/circle fcol=green size=15
d.vect end icon=basic/circle fcol=red size=15
d.vect path_1 width=1 col=blue
d.vect path_2 width=1 col=orange
d.vect path_1_new width=2 col=blue
d.vect path_2_new width=2 col=orange

Updated Shortest Path Evaluation: The thin lines correspond to the default behavior of r.walk/r.cost and r.drain, the thick lines are the updated versions by Colin.Updated Shortest Path Evaluation: The thin lines correspond to the default behavior of r.walk/r.cost and r.drain, the thick lines are the updated versions by Colin.

Dylan, your howtos are

Submitted by Stefano Costa (not verified) on Wed, 2008-02-20 09:46.

Dylan,
your howtos are always so clear and useful. Thanks for making this documentation available.