plot_rgl_model.R

require(reshape2)
require(rgl)

plot_rgl_model <- function(fdata){

  fdata <- fdata[order(fdata[, 1], fdata[, 2]), ]

  orig_names <- colnames(fdata)
  colnames(fdata) <- c("x", "y", "z")
  fdata <- as.data.frame(fdata)
  
  ## work out the min and max of x,y,z
  xlimits <- c(min(fdata$x, na.rm = T), max(fdata$x, na.rm = T))
  ylimits <- c(min(fdata$y, na.rm = T), max(fdata$y, na.rm = T))
  zlimits <- c(min(fdata$z, na.rm = T), max(fdata$z, na.rm = T))
  l <- list (x = xlimits, y = ylimits, z = zlimits)
  xyz <- do.call(expand.grid, l)
  x_boundaries <- xyz$x
  y_boundaries <- xyz$y
  z_boundaries <- xyz$z
  
  # now turn fdata into a wide format for use with the rgl.surface
  fdata[, 2] <- as.character(fdata[, 2])
  fdata[, 3] <- as.character(fdata[, 3])
  #if (verbose) print(class(fdata[, 2]))
  wide_form <- dcast(fdata, y ~ x, value.var = "z")

  wide_form_values <- as.matrix(wide_form[, 2:ncol(wide_form)])  

  x_values <- as.numeric(colnames(wide_form[2:ncol(wide_form)]))
  y_values <- as.numeric(wide_form[, 1])

  wide_form_values <- wide_form_values[order(y_values), order(x_values)]
  wide_form_values <- as.numeric(wide_form_values)
  x_values <- x_values[order(x_values)]
  y_values <- y_values[order(y_values)]
  
  zlim <- range(wide_form_values)

  zlen <- zlim[2] - zlim[1] + 1

  col <- "grey"
  
  open3d()
  plot3d(x_boundaries, y_boundaries, z_boundaries, 
         box = T, col = "black",  xlab = orig_names[1], 
         ylab = orig_names[2], zlab = orig_names[3])

  points3d(fdata, col = "blue") 

}