recexcavAAR: Vignette >>Trench visualisation<<

library(devtools)
library(recexcavAAR)
library(kriging)
library(rgl)

Like in the first vignette Ifri el Baroud we start by setting up an artificial and pretty simple excavation trench with a depth of 1.4 meters, a length of 6 meters and a width of 8 meters. This one is not parallel to the main axis of our coordinate system.

edges <- data.frame(
  x = c(6.899, 10.658, 4.428, 0.669, 6.899, 10.658, 4.428, 0.669),
  y = c(19.292, 14.616, 9.597, 14.273, 19.292, 14.616, 9.597, 14.273),
  z = c(9.7, 9.7, 9.7, 9.7, 8.3, 8.3, 8.3, 8.3)
)

We can plot the edges of this trench with rgl::plot3d, but first of all we have to calculate a reasonable aspectratio. Remember: The trench is tilted in relation to the main axis.

rangex <- abs(max(edges$x) - min(edges$x))
rangey <- abs(max(edges$y) - min(edges$y))

edgesordered = rbind(
  edges[1:4, ], 
  edges[1, ], 
  edges[5:8, ], 
  edges[5, ],
  edges[c(6,2), ],
  edges[c(3,7), ],
  edges[c(8,4), ]
)

plot3d(
  edgesordered$x, edgesordered$y, edgesordered$z,
  type="l",
  aspect = c(rangex, rangey, 5),
  xlab = "x", ylab = "y", zlab = "z",
  sub = "Grab me and rotate me!",
  col = "darkgreen"
)

bbox3d(
  xat = c(2, 4, 6, 8, 10),
  yat = c(10, 12, 14, 16, 18),
  zat = c(8.5, 9, 9.5),
  back = "lines"
)

The approach for the excavation of this trench is to go deeper in squares of 1 meter by 1 meter. The spit depth is about 30 centimeters, but as we also try to follow natural layers the individual spit depth varies a lot. Fortunately we have niveau measurements of the resulting surfaces in the dataset recexcavAAR::KT_spits.

sp <- KT_spits

splist <- list()
spitnames <- c("^surface", "^spit1", "^spit2", "^spit3", "^bottom")

for (i in 1:length(spitnames)){
  splist[[i]] <- sp[grep(spitnames[i], sp$id), ]
}

Let’s apply kriging with recexcavAAR::kriglist, transform the result with recexcavAAR::spatialwide and add the surfaces to the plot.

# I had to choose a very low pixel value to keep the vignette small enough 
maps <- kriglist(splist, x = 2, y = 3, z = 4, lags = 3, model = "spherical", pixels = 30)

surf <- list()
for (i in 1:length(maps)) {
  surf[[i]] <- spatialwide(maps[[i]]$x, maps[[i]]$y, maps[[i]]$pred, digits = 3)
}

idvec <- c()
for (i in 1:length(surf)) {
  idvec[i] <- surface3d(
    surf[[i]]$x, surf[[i]]$y, t(surf[[i]]$z),
    color = c("black", "white"),
    alpha = 0.5,
    add = TRUE
  )
}

recexcavAAR: Vignette >>Trench visualisation<<

Clemens Schmid

August 2016