# save plotting parameters
pm <- par("mfrow")
## =======================================================================
## Create topography, data
## =======================================================================
# The topographic surface
x <- seq(-pi, pi, by = 0.2)
y <- seq(0, pi, by = 0.1)
M <- mesh(x, y)
z <- with(M, sin(x)*sin(y))
# Initial condition
xi <- c(0.25 * rnorm(100) - pi/2, 0.25 * rnorm(100) - pi/4)
yi <- 0.25 * rnorm(200) + pi/2
zi <- 0.005*rnorm(200) + 0.5
# the species
species <- c(rep(1, 100), rep(2, 100))
# set initial conditions
xp <- xi; yp <- yi; zp <- zi
## =======================================================================
## Traditional graphics
## =======================================================================
par(mfrow = c(2, 2))
# Topography is defined by argument surf
for (i in 1:4) {
# update tracer distribution
xp <- xp + 0.25 * rnorm(200)
yp <- yp + 0.025 * rnorm(200)
zp <- zp + 0.25 *rnorm(200)
# plot new tracer distribution
tracers3D(xp, yp, zp, colvar = species, pch = ".", cex = 5,
main = paste("timestep ", i), col = c("orange", "blue"),
surf = list(x, y, z = z, theta = 0, facets = FALSE),
colkey = list(side = 1, length = 0.5, labels = c("sp1","sp2"),
at = c(1.25, 1.75), dist = 0.075))
}
# same, but creating topography first
if (FALSE) {
# create the topography on which to add points
persp3D(x, y, z = z, theta = 0, facets = FALSE, plot = FALSE)
for (i in 1:4) {
# update tracer distribution
xp <- xp + 0.25 * rnorm(200)
yp <- yp + 0.025 * rnorm(200)
zp <- zp + 0.25 *rnorm(200)
# plot new tracer distribution
tracers3D(xp, yp, zp, colvar = species, pch = ".", cex = 5,
main = paste("timestep ", i), col = c("orange", "blue"),
colkey = list(side = 1, length = 0.5, labels = c("sp1","sp2"),
at = c(1.25, 1.75), dist = 0.075))
}
}
## =======================================================================
## rgl graphics
## =======================================================================
# pause <- 0.05
# create a suitable topography
persp3D(x, y, z = z, theta = 0, facets = NA, plot = FALSE)
plotrgl( )
xp <- xi; yp <- yi; zp <- zi
nstep <- 10
for (i in 1:nstep) {
xp <- xp + 0.05 * rnorm(200) + 0.05
yp <- yp + 0.0025 * (rnorm(200) + 0.0025)
zp <- zp + 0.05 *rnorm(200)
# tracers3Drgl(xp, yp, zp, col = c(rep("orange", 100), rep("blue", 100)),
# main = paste("timestep ", i))
# or:
tracers3Drgl(xp, yp, zp, colvar = species, col = c("orange", "blue"),
main = paste("timestep ", i))
# Sys.sleep(pause)
# or: readline("hit enter for next")
}
# using function moviepoints3D
if (FALSE) {
# first create the data in matrices
xp <- matrix(nrow = 200, ncol = nstep, data = rep(xi, times=nstep))
yp <- matrix(nrow = 200, ncol = nstep, data = rep(yi, times=nstep))
zp <- matrix(nrow = 200, ncol = nstep, data = rep(zi, times=nstep))
tp <- matrix(nrow = 200, ncol = nstep, data = 0)
cv <- matrix(nrow = 200, ncol = nstep, data = rep(species, times=nstep))
nstep <- 10
for (i in 2:nstep) {
xp[,i] <- xp[,i-1] + 0.05 * rnorm(200) + 0.05
yp[,i] <- yp[,i-1] + 0.0025 * (rnorm(200) + 0.0025)
zp[,i] <- zp[,i-1] + 0.05 *rnorm(200)
tp[,i] <- i
}
# create the topography
persp3Drgl(x, y, z = z, theta = 0, lighting = TRUE, smooth = TRUE)
# add moviepoints:
moviepoints3D (xp, yp, zp, colvar = cv, t = tp,
wait = 0.05, cex = 10, col = c("red", "orange"))
}
# reset plotting parameters
par(mfrow = pm)
Run the code above in your browser using DataLab