Supplementary material

Spatial data analysis is a vast and rapidly growing field, and there is much that is new, and much that I did not know about when I wrote the book. As I come across new material, I will post it here. If you know of something that you feel would be good to be included, I encourage you to write and tell me about it.

Return

> thsn.sp <- as(thsn.pp, “SpatialPolygons”)

will not work. As a work-around, you can use the following functions, which are provided by Adrian Baddeley on the R-Sig-Geo list at https://stat.ethz.ch/pipermail/r-sig-geo/2009-May/005781.html.

# convert spatstat objects to sp classes
owin2Polygons <- function(x, id="1") {
  stopifnot(is.owin(x))
  x <- as.polygonal(x)
  closering <- function(df) { df[c(seq(nrow(df)), 1), ] }
  pieces <- lapply(x$bdry,
  function(p) {
  Polygon(coords=closering(cbind(p$x,p$y)),
  hole=is.hole.xypolygon(p)) })
  z <- Polygons(pieces, id)
  return(z)
}

tess2SP <- function(x) {
  stopifnot(is.tess(x))
  y <- tiles(x)
  nam <- names(y)
  z <- list()
  for(i in seq(y))
  z[[i]] <- owin2Polygons(y[[i]], nam[i])
  return(SpatialPolygons(z))
}

owin2SP <- function(x) {
  stopifnot(is.owin(x))
  y <- owin2Polygons(x)
  z <- SpatialPolygons(list(y))
  return(z)
}

Thus, for example, using these functions you would replace the second line of code on page 88 with

thsn.sp <- tess2SP(thsn.pp)

# Put legend inside map
# Get the current color settings
tr.settings <- trellis.par.get()
# Change the legend settings to the grey scale
tr.settings$superpose.polygon$col <- greys
trellis.par.set(tr.settings)
spplot(data.Set1, "HtClass2", col.regions = greys, # Fig. 7.3a
  scales = list(draw = TRUE), xlab = "Easting", ylab = "Northing",
  main = "Vegetation Height Class", sp.layout = list(obs.list),
  xlim = c(577000,578500), ylim = c(4417500,4418800), colorkey = FALSE,
  key = simpleKey(levels(data.Set1$HtClass2), x = 0.05, y = 0.05,
  corner = c(0,0), points = FALSE, rect = TRUE))

map.plot <- spplot(data.Set1, "HtClass2", col.regions = greys, # Fig. 7.3a
  scales = list(draw = TRUE), xlab = "Easting", ylab = "Northing",
  main = "Vegetation Height Class", sp.layout = list(obs.list),
  xlim = c(577000,578500), ylim = c(4417500,4418800))
names(map.plot$legend) <- "inside"
map.plot$legend$inside$x <- 0.05
map.plot$legend$inside$y <- 0.15
print(map.plot)