1

Background

I wish to calculate in R the isochrones moving from a location outwards, implementing the Tobler's hiking function via the gdistance package (https://cran.r-project.org/web/packages/gdistance/vignettes/gdistance1.pdf). An example of the implementation is provided at the above link, pp. 14-46; please note that in that reference, the issue is approached from a slighlty different perspective, i.e. least-cost path between two locations A and B.

Issue

I can't understand why I get two different results usign the same code over two dataset (the one provided by the mentioned reference vs. my toy DTM). One result make intuitively sense, the other I believe is wrong.

Reproducible code

The code below is what I used for the two analyses. The first code uses a toy DTM which comes with the gdistance package itself. The dtm and point shapefile for the second analysis are available here:

https://www.dropbox.com/s/8t0fkuovehrsf0z/mydtm.rda?dl=0

https://www.dropbox.com/s/rlzclkzlv37v2v0/loc.rda?dl=0

Following the abovementioned source (with sligth modifications), I get the first (reasonable) result attached (1st and 2nd image):

dtm <- raster(system.file("external/maungawhau.grd", package="gdistance"))
altDiff <- function(x){x[2] - x[1]}
hd <- transition(dtm, altDiff, 8, symm=FALSE)
slope <- geoCorrection(hd)
adj <- adjacent(dtm, cells=1:ncell(dtm), pairs=TRUE, directions=8)
speed <- slope
speed[adj] <- 1/(6 * exp(-3.5 * abs(slope[adj] + 0.05)))
Conductance <- geoCorrection(speed)
A <- c(2667670, 6479000)
accum <- accCost(Conductance, A)
plot(dtm)
plot(accum)
points(A[1], A[2], pch=20)

At the best of my understanding, the accum raster is (correctly, I think) representing the time it takes to move from the source location. Side note: I used the inverse of the Tobler's function to get time not speed).

Different story with my DTM (3rd and 4th image attached):

dtm <- mydtm
altDiff <- function(x){x[2] - x[1]}
hd <- transition(dtm, altDiff, 8, symm=FALSE)
slope <- geoCorrection(hd)
adj <- adjacent(dtm, cells=1:ncell(dtm), pairs=TRUE, directions=8)
speed <- slope
speed[adj] <- 1/(6 * exp(-3.5 * abs(slope[adj] + 0.05)))
Conductance <- geoCorrection(speed)
A <- cbind(loc@coords[,1], loc@coords[,2])
accum <- accCost(Conductance, A)
plot(dtm)
plot(accum)
plot(loc, add=TRUE, pch=20)

Question

At the best of my understanding, from the "exagonal-ish" outlook of raster, I have the feeling that something is going wrong in the process of accounting for the slope. I have seen something like this in case of isotropic accumulated cost-surfaces, and I do not expect to see this since the above code should implement on anisotropic cost-surface analysis (as correctly shown using the toy dtm which comes with the package). I keep getting the same "strange" result even if I re-run the analysis using a larger extent of the dtm (featuring a large variation in the elevation values). Can anyone please help me in locating what I am doing wrong?

volcano ratser DTM accum cost surf expressing time, based on vulcano raster my dtm accum cost surface expressing time, based on my dtm

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.