Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
raytracing volcano dataset in R
# Code taken from article: http://www.tylermw.com/throwing-shade/
# Author: Tyler Morgan-Wall
# License: GPLv3
volcanoshadow = matrix(1, ncol = ncol(volcano), nrow = nrow(volcano))
volc = list(x=1:nrow(volcano), y=1:ncol(volcano), z=volcano)
sunangle = 45 / 180*pi
angle = -90 / 180 * pi
diffangle = 90 / 180 * pi
numberangles = 25
anglebreaks = seq(angle, diffangle, length.out = numberangles)
maxdistance = floor(sqrt(ncol(volcano)^2 + nrow(volcano)^2))
for (i in 1:nrow(volcano)) {
for (j in 1:ncol(volcano)) {
for (anglei in anglebreaks) {
for (k in 1:maxdistance) {
xcoord = (i + sin(sunangle)*k)
ycoord = (j + cos(sunangle)*k)
if(xcoord > nrow(volcano) ||
ycoord > ncol(volcano) ||
xcoord < 0 || ycoord < 0) {
break
} else {
tanangheight = volcano[i, j] + tan(anglei) * k
if (all(c(volcano[ceiling(xcoord), ceiling(ycoord)],
volcano[floor(xcoord), ceiling(ycoord)],
volcano[ceiling(xcoord), floor(ycoord)],
volcano[floor(xcoord), floor(ycoord)]) < tanangheight)) next
if (tanangheight < akima::bilinear(volc$x,volc$y,volc$z,x0=xcoord,y0=ycoord)$z) {
volcanoshadow[i, j] = volcanoshadow[i, j] - 1 / length(anglebreaks)
break
}
}
}
}
}
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.