plot.ppp.R
#
# plot.ppp.R
#
# $Revision: 1.90 $ $Date: 2016/07/25 09:26:31 $
#
#
#--------------------------------------------------------------------------
plot.ppp <- local({
transparencyfun <- function(n) {
if(n <= 100) 1 else (0.2 + 0.8 * exp(-(n-100)/1000))
}
## determine symbol map for marks of points
default.symap.points <- function(x, ...,
chars=NULL, cols=NULL,
maxsize=NULL, meansize=NULL, markscale=NULL) {
marx <- marks(x)
if(is.null(marx)) {
## null or constant symbol map
## consider using transparent colours
if(is.null(cols) &&
!any(c("col", "fg", "bg") %in% names(list(...))) &&
(nx <- npoints(x)) > 100 &&
identical(dev.capabilities()$semiTransparency, TRUE) &&
spatstat.options("transparent"))
cols <- rgb(0,0,0,transparencyfun(nx))
return(symbolmap(..., chars=chars, cols=cols))
}
if(!is.null(dim(marx)))
stop("Internal error: multivariate marks in default.symap.points")
argnames <- names(list(...))
shapegiven <- "shape" %in% argnames
chargiven <- (!is.null(chars)) || ("pch" %in% argnames)
assumecircles <- !(shapegiven || chargiven)
sizegiven <- ("size" %in% argnames) ||
(("cex" %in% argnames) && !shapegiven)
if(inherits(marx, c("Date", "POSIXt"))) {
## ......... marks are dates or date/times .....................
timerange <- range(marx, na.rm=TRUE)
shapedefault <- if(!assumecircles) list() else list(shape="circles")
if(sizegiven) {
g <- do.call(symbolmap,
resolve.defaults(list(range=timerange),
list(...),
shapedefault,
list(chars=chars, cols=cols)))
return(g)
}
## attempt to determine a scale for the marks
y <- scaletointerval(marx, 0, 1, timerange)
y <- y[is.finite(y)]
if(length(y) == 0) return(symbolmap(..., chars=chars, cols=cols))
scal <- mark.scale.default(y, as.owin(x),
markscale=markscale, maxsize=maxsize,
meansize=meansize,
characters=chargiven)
if(is.na(scal)) return(symbolmap(..., chars=chars, cols=cols))
## scale determined
sizefun <- function(x, scal=1) {
(scal/2) * scaletointerval(x, 0, 1, timerange)
}
formals(sizefun)[[2]] <- scal ## ensures value of 'scal' is printed
##
g <- do.call(symbolmap,
resolve.defaults(list(range=timerange),
list(...),
shapedefault,
list(size=sizefun)))
return(g)
}
if(is.numeric(marx)) {
## ............. marks are numeric values ...................
marx <- marx[is.finite(marx)]
if(length(marx) == 0)
return(symbolmap(..., chars=chars, cols=cols))
markrange <- range(marx)
##
if(sizegiven) {
g <- do.call(symbolmap,
resolve.defaults(list(range=markrange),
list(...),
if(assumecircles) list(shape="circles") else list(),
list(chars=chars, cols=cols)))
return(g)
}
## attempt to determine a scale for the marks
if(all(markrange == 0))
return(symbolmap(..., chars=chars, cols=cols))
scal <- mark.scale.default(marx, as.owin(x),
markscale=markscale, maxsize=maxsize,
meansize=meansize,
characters=chargiven)
if(is.na(scal)) return(symbolmap(..., chars=chars, cols=cols))
## scale determined
if(markrange[1] >= 0) {
## all marks are nonnegative
shapedefault <-
if(!assumecircles) list() else list(shape="circles")
cexfun <- function(x, scal=1) { scal * x }
circfun <- function(x, scal=1) { scal * x }
formals(cexfun)[[2]] <- formals(circfun)[[2]] <- scal
sizedefault <-
if(sizegiven) list() else
if(chargiven) list(cex=cexfun) else list(size=circfun)
} else {
## some marks are negative
shapedefault <-
if(!assumecircles) list() else
list(shape=function(x) { ifelse(x >= 0, "circles", "squares") })
cexfun <- function(x, scal=1) { scal * abs(x) }
circfun <- function(x, scal=1) { scal * abs(x) }
formals(cexfun)[[2]] <- formals(circfun)[[2]] <- scal
sizedefault <-
if(sizegiven) list() else
if(chargiven) list(cex=cexfun) else list(size=circfun)
}
g <- do.call(symbolmap,
resolve.defaults(list(range=markrange),
list(...),
shapedefault,
sizedefault,
list(chars=chars, cols=cols)))
return(g)
}
## ........... non-numeric marks .........................
um <- if(is.factor(marx)) levels(marx) else sort(unique(marx))
ntypes <- length(um)
if(!is.null(cols))
cols <- rep.int(cols, ntypes)[1:ntypes]
if(shapegiven && sizegiven) {
#' values mapped to symbols (shape and size specified)
g <- symbolmap(inputs=um, ..., cols=cols)
} else if(!shapegiven) {
#' values mapped to 'pch'
chars <- default.charmap(ntypes, chars)
g <- symbolmap(inputs=um, ..., chars=chars, cols=cols)
} else {
#' values mapped to symbols
#' determine size
scal <- mark.scale.default(rep(1, npoints(x)),
Window(x),
maxsize=maxsize,
meansize=meansize,
characters=FALSE)
g <- symbolmap(inputs=um, ..., size=scal, cols=cols)
}
return(g)
}
default.charmap <- function(n, ch=NULL) {
if(!is.null(ch))
return(rep.int(ch, n)[1:n])
if(n <= 25)
return(1:n)
ltr <- c(letters, LETTERS)
if(n <= 52)
return(ltr[1:n])
## wrapped sequence of letters
warning("Too many types to display every type as a different character")
return(ltr[1 + (0:(n - 1) %% 52)])
}
## main function
plot.ppp <-
function(x, main, ..., clipwin=NULL,
chars=NULL, cols=NULL, use.marks=TRUE,
which.marks=NULL, add=FALSE, type=c("p", "n"),
legend=TRUE, leg.side=c("left", "bottom", "top", "right"),
leg.args=list(),
symap=NULL, maxsize=NULL, meansize=NULL, markscale=NULL, zap=0.01,
show.window=show.all, show.all=!add, do.plot=TRUE,
multiplot=TRUE)
{
if(missing(main))
main <- short.deparse(substitute(x))
type <- match.arg(type)
if(!missing(maxsize) || !missing(markscale) || !missing(meansize))
warn.once("circlescale",
"Interpretation of arguments maxsize and markscale",
"has changed (in spatstat version 1.37-0 and later).",
"Size of a circle is now measured by its diameter.")
if(clipped <- !is.null(clipwin)) {
stopifnot(is.owin(clipwin))
W <- Window(x)
clippy <- if(is.mask(W)) intersect.owin(W, clipwin) else edges(W)[clipwin]
x <- x[clipwin]
} else clippy <- NULL
## sensible default position
legend <- legend && show.all
if(legend) {
leg.side <- match.arg(leg.side)
vertical <- (leg.side %in% c("left", "right"))
}
# if(type == "n" || npoints(x) == 0) {
# ## plot the window only
# xwindow <- x$window
# if(do.plot)
# do.call(plot.owin,
# resolve.defaults(list(xwindow),
# list(...),
# list(main=main, invert=TRUE, add=add,
# type=if(show.window) "w" else "n")))
# if(is.null(symap)) symap <- symbolmap()
# attr(symap, "bbox") <- as.rectangle(xwindow)
# return(invisible(symap))
# }
## ................................................................
## Handle multiple columns of marks as separate plots
## (unless add=TRUE or which.marks selects a single column
## or multipage = FALSE)
if(use.marks && is.data.frame(mx <- marks(x))) {
implied.all <- is.null(which.marks)
want.several <- implied.all || is.data.frame(mx <- mx[,which.marks])
do.several <- want.several && !add && multiplot
if(do.several) {
## generate one plot for each column of marks
y <- solapply(mx, setmarks, x=x)
out <- do.call(plot,
resolve.defaults(list(x=y, main=main,
show.window=show.window && !clipped,
do.plot=do.plot,
type=type),
list(...),
list(equal.scales=TRUE),
list(panel.end=clippy),
list(legend=legend,
leg.side=leg.side,
leg.args=leg.args),
list(chars=chars, cols=cols,
maxsize=maxsize,
meansize=meansize,
markscale=markscale,
zap=zap)))
return(invisible(out))
}
if(is.null(which.marks)) {
which.marks <- 1
if(do.plot) message("Plotting the first column of marks")
}
}
## ............... unmarked, or single column of marks ....................
## Determine symbol map and mark values to be used
y <- x
if(!is.marked(x, na.action="ignore") || !use.marks) {
## Marks are not mapped.
marx <- NULL
if(is.null(symap))
symap <- default.symap.points(unmark(x), ..., chars=chars, cols=cols)
} else {
## Marked point pattern
marx <- marks(y, dfok=TRUE)
if(is.data.frame(marx)) {
## select column or take first colum
marx <- marx[, which.marks]
y <- setmarks(y, marx)
}
if(npoints(y) > 0) {
ok <- complete.cases(as.data.frame(y))
if(!any(ok)) {
warning("All mark values are NA; plotting locations only.")
if(is.null(symap))
symap <- default.symap.points(unmark(x), ..., chars=chars, cols=cols)
} else if(any(!ok)) {
warning(paste("Some marks are NA;",
"corresponding points are omitted."))
x <- x[ok]
y <- y[ok]
marx <- marks(y)
}
}
## apply default symbol map
if(is.null(symap))
symap <- default.symap.points(y, chars=chars, cols=cols,
maxsize=maxsize, meansize=meansize,
markscale=markscale,
...)
}
# gtype <- symbolmaptype(symap)
## Determine bounding box for main plot
BB <- as.rectangle(x)
sick <- inherits(x, "ppp") && !is.null(rejects <- attr(x, "rejects"))
if(sick) {
## Get relevant parameters
par.direct <- list(main=main, use.marks=use.marks,
maxsize=maxsize, meansize=meansize, markscale=markscale)
par.rejects <- resolve.1.default(list(par.rejects=list(pch="+")),
list(...))
par.all <- resolve.defaults(par.rejects, par.direct)
rw <- resolve.defaults(list(...), list(rejectwindow=NULL))$rejectwindow
## determine window for rejects
rwin <-
if(is.null(rw))
rejects$window
else if(is.logical(rw) && rw)
rejects$window
else if(inherits(rw, "owin"))
rw
else if(is.character(rw)) {
switch(rw,
box={boundingbox(rejects, x)},
ripras={ripras(c(rejects$x, x$x), c(rejects$y, x$y))},
stop(paste("Unrecognised option: rejectwindow=", rw)))
} else stop("Unrecognised format for rejectwindow")
if(is.null(rwin))
stop("Selected window for rejects pattern is NULL")
BB <- boundingbox(BB, as.rectangle(rwin))
}
## Augment bounding box with space for legend, if appropriate
legend <- legend && (symbolmaptype(symap) != "constant")
if(legend) {
## guess maximum size of symbols
maxsize <- invoke.symbolmap(symap, marx,
corners(as.rectangle(x)),
add=add, do.plot=FALSE)
sizeguess <- if(maxsize <= 0) NULL else (1.5 * maxsize)
leg.args <- append(list(side=leg.side, vertical=vertical), leg.args)
## draw up layout
legbox <- do.call.matched(plan.legend.layout,
append(list(B=BB, size = sizeguess,
started=FALSE, map=symap),
leg.args))
## bounding box for everything
BB <- legbox$A
}
## return now if not plotting
attr(symap, "bbox") <- BB
if(!do.plot)
return(invisible(symap))
## ............. start plotting .......................
pt <- prepareTitle(main)
main <- pt$main
nlines <- pt$nlines
blankmain <- if(nlines == 0) "" else rep(" ", nlines)
rez <- resolve.defaults(list(...),
list(cex.main=1,
xlim=NULL,
ylim=NULL))
plot(BB, type="n", add=add, main=blankmain, show.all=show.all,
cex.main=rez$cex.main, xlim=rez$xlim, ylim=rez$ylim)
if(sick) {
if(show.window) {
## plot windows
if(!is.null(rw)) {
## plot window for rejects
rwinpardefault <- list(lty=2,lwd=1,border=1)
rwinpars <-
resolve.defaults(par.rejects, rwinpardefault)[names(rwinpardefault)]
do.call(plot.owin, append(list(rwin, add=TRUE), rwinpars))
}
## plot window of main pattern
if(!clipped) {
do.call(plot.owin,
resolve.defaults(list(x$window, add=TRUE),
list(...),
list(invert=TRUE)))
} else plot(clippy, add=TRUE, ...)
}
if(type != "n") {
## plot reject points
do.call(plot.ppp, append(list(rejects, add=TRUE), par.all))
warning(paste(rejects$n, "illegal points also plotted"))
}
## the rest is added
add <- TRUE
}
## Now convert to bona fide point pattern
x <- as.ppp(x)
xwindow <- x$window
## Plot observation window (or at least the main title)
do.call(plot.owin,
resolve.defaults(list(x=xwindow,
add=TRUE,
main=main,
type=if(show.window && !clipped) "w" else "n",
show.all=show.all),
list(...),
list(invert=TRUE)))
## If clipped, plot visible part of original window
if(show.window && clipped)
plot(clippy, add=TRUE, ...)
# else if(show.all) fakemaintitle(as.rectangle(xwindow), main, ...)
if(type != "n") {
## plot symbols ##
invoke.symbolmap(symap, marx, x, add=TRUE)
}
## add legend
if(legend) {
b <- legbox$b
legendmap <- if(length(leg.args) == 0) symap else
do.call(update, append(list(object=symap), leg.args))
do.call(plot,
append(list(x=legendmap, main="", add=TRUE,
xlim=b$xrange, ylim=b$yrange),
leg.args))
}
return(invisible(symap))
}
plot.ppp
})
mark.scale.default <- function(marx, w, markscale=NULL,
maxsize=NULL, meansize=NULL,
characters=FALSE) {
## establish values of markscale, maxsize, meansize
ngiven <- (!is.null(markscale)) +
(!is.null(maxsize)) +
(!is.null(meansize))
if(ngiven > 1)
stop("Only one of the arguments markscale, maxsize, meansize",
" should be given", call.=FALSE)
if(ngiven == 0) {
## if ALL are absent, enforce the spatstat defaults
## (which could also be null)
pop <- spatstat.options("par.points")
markscale <- pop$markscale
maxsize <- pop$maxsize
meansize <- pop$meansize
}
## Now check whether markscale is fixed
if(!is.null(markscale)) {
stopifnot(markscale > 0)
return(markscale)
}
# Usual case: markscale is to be determined from maximum/mean physical size
if(is.null(maxsize) && is.null(meansize)) {
## compute default value of 'maxsize'
## guess appropriate max physical size of symbols
bb <- as.rectangle(w)
maxsize <- 1.4/sqrt(pi * length(marx)/area(bb))
maxsize <- min(maxsize, diameter(bb) * 0.07)
## updated: maxsize now represents *diameter*
maxsize <- 2 * maxsize
} else {
if(!is.null(maxsize)) stopifnot(maxsize > 0) else stopifnot(meansize > 0)
}
# Examine mark values
absmarx <- abs(marx)
maxabs <- max(absmarx)
tiny <- (maxabs < 4 * .Machine$double.eps)
if(tiny)
return(NA)
## finally determine physical scale for symbols
if(!is.null(maxsize)) {
scal <- maxsize/maxabs
} else {
meanabs <- mean(absmarx)
scal <- meansize/meanabs
}
if(!characters) return(scal)
## if using characters ('pch') we need to
## convert physical sizes to 'cex' values
charsize <- max(sidelengths(as.rectangle(w)))/40
return(scal/charsize)
}
fakemaintitle <- function(bb, main, ...) {
## Try to imitate effect of 'title(main=main)' above a specified box
if(!any(nzchar(main))) return(invisible(NULL))
bb <- as.rectangle(bb)
x0 <- mean(bb$xrange)
y0 <- bb$yrange[2] + length(main) * diff(bb$yrange)/12
parnames <- c('cex.main', 'col.main', 'font.main')
parlist <- par(parnames)
parlist <- resolve.defaults(list(...), parlist)[parnames]
names(parlist) <- c('cex', 'col', 'font')
do.call.matched(text.default,
resolve.defaults(list(x=x0, y=y0, labels=main),
parlist, list(...)),
funargs=graphicsPars("text"))
return(invisible(NULL))
}