Nothing
#
# tess.R
#
# support for tessellations
#
# $Revision: 1.111 $ $Date: 2024/02/04 08:04:51 $
#
tess <- function(..., xgrid=NULL, ygrid=NULL, tiles=NULL, image=NULL,
window=NULL, marks=NULL, keepempty=FALSE,
unitname=NULL, check=TRUE) {
uname <- unitname
if(!is.null(window)) {
window <- as.owin(window)
if(is.null(uname)) uname <- unitname(window)
}
isrect <- !is.null(xgrid) && !is.null(ygrid)
istiled <- !is.null(tiles)
isimage <- !is.null(image)
if(isrect + istiled + isimage != 1)
stop("Must specify either (xgrid, ygrid) or tiles or img")
if(isrect) {
stopifnot(is.numeric(xgrid) && all(diff(xgrid) > 0))
stopifnot(is.numeric(ygrid) && all(diff(ygrid) > 0))
if(!is.null(window))
warning("Argument 'window' ignored, because xgrid, grid are given")
window <- owinInternalRect(range(xgrid), range(ygrid), unitname=uname)
ntiles <- (length(xgrid)-1) * (length(ygrid)-1)
out <- list(type="rect", window=window, xgrid=xgrid, ygrid=ygrid, n=ntiles)
} else if(istiled) {
stopifnot(is.list(tiles))
if(check) {
if(!all(sapply(tiles, is.owin)))
stop("Tiles must be a list of owin objects")
if(!is.null(uname)) {
## attach new unit name to each tile
tiles <- lapply(tiles, "unitname<-", value=uname)
} else {
## extract unit names from tiles, check agreement, use as unitname
uu <- unique(lapply(tiles, unitname))
uu <- uu[!sapply(uu, is.null)]
nun <- length(uu)
if(nun > 1)
stop("Tiles have inconsistent names for the unit of length")
if(nun == 1) {
## use this unit name
uname <- uu[[1]]
if(!is.null(window))
unitname(window) <- uname
}
}
}
if(!keepempty && check) {
# remove empty tiles
isempty <- sapply(tiles, is.empty)
if(all(isempty))
stop("All tiles are empty")
if(any(isempty))
tiles <- tiles[!isempty]
}
ntiles <- length(tiles)
nam <- names(tiles)
lev <- if(!is.null(nam) && all(nzchar(nam))) nam else 1:ntiles
if(is.null(window))
window <- do.call(union.owin, unname(tiles))
if(is.mask(window) || any(sapply(tiles, is.mask))) {
# convert to pixel image tessellation
Grid <- do.call(commonGrid, append(list(window), unname(tiles)))
ima <- as.im(window, W=Grid)
ima$v[] <- NA
for(i in 1:ntiles)
ima[tiles[[i]]] <- i
ima <- ima[window, drop=FALSE]
ima <- eval.im(factor(ima, levels=1:ntiles))
levels(ima) <- lev
out <- list(type="image",
window=window, image=ima, n=length(lev))
} else {
# tile list
window <- rescue.rectangle(window)
out <- list(type="tiled", window=window, tiles=tiles, n=length(tiles))
}
} else if(isimage) {
# convert to factor valued image
image <- as.im(image)
if(!is.null(uname)) unitname(image) <- uname
switch(image$type,
logical={
# convert to factor
if(keepempty)
image <- eval.im(factor(image, levels=c(FALSE,TRUE)))
else
image <- eval.im(factor(image))
},
factor={
# eradicate unused levels
if(!keepempty)
image <- eval.im(factor(image))
},
{
# convert to factor
image <- eval.im(factor(image))
})
if(is.null(window)) window <- as.owin(image)
out <- list(type="image", window=window, image=image, n=length(levels(image)))
} else stop("Internal error: unrecognised format")
## add marks!
if(!is.null(marks)) {
marks <- as.data.frame(marks)
if(nrow(marks) != out$n)
stop(paste("wrong number of marks:",
nrow(marks), "should be", out$n),
call.=FALSE)
out$marks <- marks
}
class(out) <- c("tess", class(out))
return(out)
}
is.tess <- function(x) { inherits(x, "tess") }
print.tess <- function(x, ..., brief=FALSE) {
full <- !brief
if(full) cat("Tessellation\n")
win <- x$window
switch(x$type,
rect={
if(full) {
unitinfo <- summary(unitname(win))
if(evenly.spaced(x$xgrid) && evenly.spaced(x$ygrid))
splat("Tiles are equal rectangles, of dimension",
signif(mean(diff(x$xgrid)), 5),
"x",
signif(mean(diff(x$ygrid)), 5),
unitinfo$plural, " ", unitinfo$explain)
else
splat("Tiles are unequal rectangles")
}
splat(length(x$xgrid)-1, "by", length(x$ygrid)-1, "grid of tiles")
},
tiled={
if(full) {
if(win$type == "polygonal")
splat("Tiles are irregular polygons")
else
splat("Tiles are windows of general type")
}
splat(length(x$tiles), "tiles (irregular windows)")
},
image={
nlev <- length(levels(x$image))
if(full) {
splat("Tessellation is determined by a factor-valued image with",
nlev, "levels")
} else splat(nlev, "tiles (levels of a pixel image)")
})
if(!is.null(marx <- x$marks)) {
m <- dim(marx)[2] %orifnull% 1
if(m == 1) splat("Tessellation is marked") else
splat("Tessellation has", m, "columns of marks:",
commasep(sQuote(colnames(marx))))
}
if(full) print(win)
invisible(NULL)
}
unitname.tess <- function(x) unitname(x$window)
"unitname<-.tess" <- function(x, value) {
unitname(x$window) <- value
switch(x$type,
rect={},
tiled={
x$tiles <- lapply(x$tiles, "unitname<-", value)
},
image={
unitname(x$image) <- value
})
return(x)
}
plot.tess <- local({
plotpars <- c("sub", "lty", "lwd",
"cex.main", "col.main", "font.main",
"cex.sub", "col.sub", "font.sub", "border")
plot.tess <- function(x, ..., main, add=FALSE, show.all=!add,
border=NULL,
do.plot=TRUE,
do.labels=FALSE, labels=tilenames(x),
labelargs=list(),
do.col=FALSE,
values=marks(x),
multiplot=TRUE,
col=NULL,
ribargs=list()) {
if(missing(main) || is.null(main))
main <- short.deparse(substitute(x))
ntiles <- x$n
if(!do.col) {
#' Plot tiles, with adornment
y <- NULL
result <- NULL
bbox <- NULL
need.legend <- FALSE
} else {
#' Fill tiles with colours determined by 'values'
if(markformat(values) == "hyperframe")
values <- as.data.frame(values) #' automatic warning
#' Determine values associated with each tile
switch(markformat(values),
none = {
#' no values assigned.
#' default is tile name
tn <- tilenames(x)
values <- factor(tn, levels=tn)
},
vector = {
#' vector of values.
#' validate length of vector
check.anyvector(values, ntiles, things="tiles", vname="values")
},
dataframe = {
#' data frame or matrix of values.
values <- as.data.frame(values)
if(nrow(values) != ntiles)
stop(paste("Number of rows of values =", nrow(values),
"!=", ntiles, "= number of tiles"),
call.=FALSE)
if(multiplot && ncol(values) > 1 && !add) {
#' Multiple Panel Plot
result <- multi.plot.tess(x, ...,
main=main, show.all=show.all,
border=border, do.plot=do.plot,
do.labels=do.labels, labels=labels,
labelargs=labelargs, do.col=do.col,
col=col, ribargs=ribargs)
return(invisible(result))
}
if(ncol(values) > 1)
warning("Using only the first column of values")
values <- values[,1]
},
stop("Format of values is not understood")
)
#' Single Panel Plot
#' Determine colour map and plan layout (including colour ribbon)
#' using rules for pixel images
y <- as.im(as.function(x, values=values))
dont.complain.about(y)
result <- do.call(plot.im,
resolve.defaults(
list(x=quote(y),
do.plot=FALSE,
show.all=show.all, add=add, main=main,
col=col, ribargs=ribargs),
list(...),
list(valuesAreColours=FALSE)
))
#' exit if not actually plotting
if(!do.plot) return(invisible(result))
#' extract info
colmap <- result
bbox <- attr(result, "bbox")
bbox.legend <- attr(result, "bbox.legend")
need.legend <- !is.null(bbox.legend)
}
#' Start Plot
#' initialise plot region if it is determined
if(do.plot && !is.null(bbox) && !add) {
plot(bbox, main=" ", type="n")
add <- TRUE
}
switch(x$type,
rect={
win <- x$window
dont.complain.about(win)
z <- do.call.matched(plot.owin,
resolve.defaults(list(x=quote(win),
main=main,
add=add,
show.all=show.all,
do.plot=do.plot),
list(...)),
extrargs=plotpars)
if(is.null(result))
result <- z
if(do.plot) {
#' actually plot
if(do.col) {
#' fill tiles with colours
colours <- colmap(values)
til <- tiles(x)
for(i in seq_len(x$n))
plot(til[[i]], add=TRUE,
col=colours[i], border=border,
main="", ...)
} else {
#' draw tile boundaries only
xg <- x$xgrid
yg <- x$ygrid
do.call.matched(segments,
resolve.defaults(list(x0=xg, y0=win$yrange[1],
x1=xg, y1=win$yrange[2]),
list(col=border),
list(...),
.StripNull=TRUE))
do.call.matched(segments,
resolve.defaults(list(x0=win$xrange[1], y0=yg,
x1=win$xrange[2], y1=yg),
list(col=border),
list(...),
.StripNull=TRUE))
}
}
},
tiled={
xwin <- x$window
dont.complain.about(xwin)
z <- do.call.matched(plot.owin,
resolve.defaults(list(x=quote(xwin),
main=main,
add=add,
show.all=show.all,
do.plot=do.plot),
list(...)),
extrargs=plotpars)
if(is.null(result)) result <- z
if(do.plot) {
#' plot each tile
til <- tiles(x)
if(!do.col) {
#' border only
lapply(til, plot.owin, ..., add=TRUE, border=border)
} else {
#' fill with colour
colours <- colmap(values)
mapply(plot.owin, x=til, col=colours,
MoreArgs=list(add=TRUE, main="", border=border, ...))
}
}
},
image={
if(is.null(y)) y <- x$image
dont.complain.about(y)
result <-
do.call(plot,
resolve.defaults(list(quote(y), add=add, main=main,
show.all=show.all,
do.plot=do.plot,
col=col, ribargs=ribargs),
list(...),
list(valuesAreColours=FALSE)))
need.legend <- FALSE
})
if(do.plot && do.labels) {
labels <- paste(as.vector(labels))
til <- tiles(x)
incircles <- lapply(til, incircle)
x0 <- sapply(incircles, getElement, name="x")
y0 <- sapply(incircles, getElement, name="y")
do.call.matched(text.default,
resolve.defaults(list(x=x0, y = y0),
list(labels=labels),
labelargs),
funargs=graphicsPars("text"))
}
if(do.plot && need.legend) {
#' determine position of legend
xlim <- bbox.legend$xrange
ylim <- bbox.legend$yrange
sidecode <- attr(colmap, "side.legend")
vertical <- sidecode %in% c(2,4)
do.call(plot.colourmap,
resolve.defaults(list(x=quote(colmap),
add=TRUE, main="",
xlim=xlim, ylim=ylim,
side=sidecode, vertical=vertical),
ribargs,
list(...)))
}
return(invisible(result))
}
multi.plot.tess <- function(x, ..., zlim=NULL, col=NULL, equal.ribbon=FALSE) {
if(equal.ribbon && is.null(zlim) && !inherits(col, "colourmap"))
zlim <- range(marks(x))
if(!is.null(zlim)) {
result <- plot(unstack(x), ..., zlim=zlim, col=col)
} else {
result <- plot(unstack(x), ..., col=col)
}
return(invisible(result))
}
plot.tess
})
"[<-.tess" <- function(x, i, ..., value) {
switch(x$type,
rect=,
tiled={
til <- tiles(x)
til[i] <- value
ok <- !unlist(lapply(til, is.null))
x <- tess(tiles=til[ok])
},
image={
stop("Cannot assign new values to subsets of a pixel image")
})
return(x)
}
"[.tess" <- function(x, i, ...) {
trap.extra.arguments(..., .Context="in [.tess")
if(missing(i)) return(x)
if(is.owin(i))
return(intersect.tess(x, i))
switch(x$type,
rect=,
tiled={
til <- tiles(x)[i]
return(tess(tiles=til))
},
image={
img <- x$image
oldlev <- levels(img)
newlev <- unique(oldlev[i])
img <- eval.im(factor(img, levels=newlev))
return(tess(image=img))
})
}
tiles <- function(x) {
stopifnot(is.tess(x))
switch(x$type,
rect={
out <- list()
xg <- x$xgrid
yg <- x$ygrid
nx <- length(xg) - 1
ny <- length(yg) - 1
for(j in rev(seq_len(ny))) {
for(i in seq_len(nx)) {
winij <- owinInternalRect(xg[c(i,i+1)], yg[c(j,j+1)])
out <- append(out, list(winij))
}
}
},
tiled={
out <- x$tiles
},
image={
out <- list()
ima <- x$image
lev <- levels(ima)
for(i in seq_along(lev))
out[[i]] <- solutionset(ima == lev[i])
})
names(out) <- tilenames(x)
out <- as.solist(out)
return(out)
}
tiles.empty <- function(x) {
stopifnot(is.tess(x))
switch(x$type,
rect = {
nx <- length(x$xgrid) - 1
ny <- length(x$ygrid) - 1
ans <- rep(FALSE, nx * ny)
},
tiled = {
ans <- sapply(x$tiles, is.empty)
},
image = {
ans <- (table(x$image[]) == 0)
})
return(ans)
}
tilenames <- function(x) {
UseMethod("tilenames")
}
tilenames.tess <- function(x) {
switch(x$type,
rect={
if(!is.null(x$tilenames)) {
out <- x$tilenames
} else {
nx <- length(x$xgrid) - 1
ny <- length(x$ygrid) - 1
ij <- expand.grid(1:nx, 1:ny)
out <- paste0("Tile row ", ij[,2], ", col ", ij[,1])
}
},
tiled={
out <- names(x$tiles)
if(sum(nzchar(out)) != x$n)
out <- paste("Tile", seq_len(x$n))
},
image={
out <- levels(x$image)
}
)
return(as.character(out))
}
"tilenames<-" <- function(x, value) {
UseMethod("tilenames<-")
}
"tilenames<-.tess" <- function(x, value) {
if(!is.null(value)) {
## validate length
value <- as.character(value)
nv <- length(value)
switch(x$type,
rect = {
nx <- length(x$xgrid) - 1
ny <- length(x$ygrid) - 1
n <- nx * ny
},
tiled = { n <- length(x$tiles) },
image = { n <- length(levels(x$image)) })
if(nv != n)
stop("Replacement value has wrong length",
paren(paste(nv, "instead of", n)))
}
switch(x$type,
rect={
x$tilenames <- value
},
tiled={
names(x$tiles) <- value
},
image={
levels(x$image) <- value %orifnull% (1:n)
}
)
return(x)
}
marks.tess <- function(x, ...) {
stopifnot(is.tess(x))
return(x$marks)
}
"marks<-.tess" <- function(x, ..., value) {
stopifnot(is.tess(x))
if(!is.null(value)) {
value <- as.data.frame(value)
ntil <- x$n
if(nrow(value) != ntil)
stop(paste("replacement value for marks has wrong length:",
nrow(value), "should be", ntil),
call.=FALSE)
rownames(value) <- NULL
if(ncol(value) == 1) colnames(value) <- "marks"
}
x$marks <- value
return(x)
}
unmark.tess <- function(X) { marks(X) <- NULL; return(X) }
tile.areas <- function(x) {
stopifnot(is.tess(x))
switch(x$type,
rect={
xg <- x$xgrid
yg <- x$ygrid
a <- outer(rev(diff(yg)), diff(xg), "*")
a <- as.vector(t(a))
names(a) <- as.vector(t(tilenames(x)))
},
tiled={
a <- as.numeric(sapply(x$tiles, area))
names(a) <- tilenames(x)
},
image={
z <- x$image
a <- as.numeric(table(z$v)) * z$xstep * z$ystep
names(a) <- tilenames(x)
})
return(a)
}
as.im.tess <- function(X, W=NULL, ...,
eps=NULL, dimyx=NULL, xy=NULL,
rule.eps=c("adjust.eps", "grow.frame", "shrink.frame"),
na.replace=NULL, values=NULL) {
rule.eps <- match.arg(rule.eps)
## if W is present, it may have to be converted
if(!is.null(W)) {
stopifnot(is.owin(W))
if(W$type != "mask")
W <- as.mask(W, eps=eps, dimyx=dimyx, xy=xy, rule.eps=rule.eps)
}
if(!is.null(values)) {
if(!is.atomic(values))
stop("Argument 'values' should contain numeric, logical or factor values",
call.=FALSE)
if(length(values) != nobjects(X))
stop(paste("length(values) =", length(values), "!=",
nobjects(X), "= number of tiles"),
call.=FALSE)
}
switch(X$type,
image={
if(is.null(values)) {
## result is factor image
out <- as.im(X$image, W=W, eps=eps, dimyx=dimyx, xy=xy,
rule.eps=rule.eps,
na.replace=na.replace)
} else {
## map tiles to 'values'
out <- as.im(X$image, W=W, eps=eps, dimyx=dimyx, xy=xy,
rule.eps=rule.eps)
out <- eval.im(values[as.integer(out)])
## replace NA by other value
if(!is.null(na.replace))
out <- as.im(out, na.replace=na.replace)
}
},
tiled={
if(is.null(W))
W <- as.mask(as.owin(X), eps=eps, dimyx=dimyx, xy=xy,
rule.eps=rule.eps)
til <- X$tiles
ntil <- length(til)
nama <- names(til)
if(is.null(nama) || !all(nzchar(nama)))
nama <- paste(seq_len(ntil))
xy <- list(x=W$xcol, y=W$yrow)
for(i in seq_len(ntil)) {
indic <- as.mask(til[[i]], xy=xy)
tag <- as.im(indic, value=i)
if(i == 1) {
out <- tag
outv <- as.integer(out$v)
} else {
outv <- pmin.int(outv, tag$v, na.rm=TRUE)
}
}
if(is.null(values)) {
## result is factor image
outv <- factor(outv, levels=seq_len(ntil), labels=nama)
} else {
## map tiles to 'values'
outv <- values[outv]
}
if(!is.null(na.replace) && anyNA(outv))
outv[is.na(outv)] <- na.replace
out <- im(outv, out$xcol, out$yrow, unitname=unitname(W))
},
rect={
out <- as.im(W %orifnull% as.rectangle(X),
eps=eps, dimyx=dimyx, xy=xy,
rule.eps=rule.eps)
xg <- X$xgrid
yg <- X$ygrid
nrows <- length(yg) - 1L
ncols <- length(xg) - 1L
jx <- findInterval(out$xcol, xg, rightmost.closed=TRUE)
iy <- findInterval(out$yrow, yg, rightmost.closed=TRUE)
M <- as.matrix(out)
Jcol <- jx[col(M)]
Irow <- nrows - iy[row(M)] + 1L
Ktile <- as.integer(Jcol + ncols * (Irow - 1L))
if(is.null(values)) {
## result is factor image
outv <- factor(Ktile, levels=seq_len(nrows * ncols))
} else {
## map tiles to 'values'
outv <- values[Ktile]
}
if(!is.null(na.replace) && anyNA(outv))
outv[is.na(outv)] <- na.replace
out <- im(outv, xcol=out$xcol, yrow=out$yrow,
unitname=unitname(W))
}
)
return(out)
}
nobjects.tess <- function(x) {
switch(x$type,
image = length(levels(x$image)),
rect = (length(x$xgrid)-1) * (length(x$ygrid)-1),
tiled = length(x$tiles))
}
as.function.tess <- function(x, ..., values=NULL) {
V <- x
if(is.null(values)) {
f <- function(x,y) { tileindex(x,y,V) }
} else {
if(is.data.frame(values)) values <- unlist(values)
if(length(values) != nobjects(x))
stop("Length of 'values' should equal the number of tiles", call.=FALSE)
f <- function(x,y) { values[as.integer(tileindex(x,y,V))] }
}
g <- funxy(f, Window(V))
return(g)
}
tileindex <- function(x, y, Z) {
stopifnot(is.tess(Z))
if((missing(y) || is.null(y)) && all(c("x", "y") %in% names(x))) {
y <- x$y
x <- x$x
}
stopifnot(length(x) == length(y))
switch(Z$type,
rect={
jx <- findInterval(x, Z$xgrid, rightmost.closed=TRUE)
iy <- findInterval(y, Z$ygrid, rightmost.closed=TRUE)
nrows <- length(Z$ygrid) - 1
ncols <- length(Z$xgrid) - 1
iy[iy < 1 | iy > nrows] <- NA
jx[jx < 1 | jx > ncols] <- NA
jcol <- jx
irow <- nrows - iy + 1
ktile <- jcol + ncols * (irow - 1)
m <- factor(ktile, levels=seq_len(nrows*ncols))
ij <- expand.grid(j=seq_len(ncols),i=seq_len(nrows))
levels(m) <- paste("Tile row ", ij$i, ", col ", ij$j, sep="")
},
tiled={
n <- length(x)
todo <- seq_len(n)
nt <- length(Z$tiles)
m <- integer(n)
for(i in 1:nt) {
ti <- Z$tiles[[i]]
hit <- inside.owin(x[todo], y[todo], ti)
if(any(hit)) {
m[todo[hit]] <- i
todo <- todo[!hit]
}
if(length(todo) == 0)
break
}
m[m == 0] <- NA
nama <- names(Z$tiles)
lev <- seq_len(nt)
lab <- if(!is.null(nama) && all(nzchar(nama))) nama else paste("Tile", lev)
m <- factor(m, levels=lev, labels=lab)
},
image={
Zim <- Z$image
m <- lookup.im(Zim, x, y, naok=TRUE)
if(anyNA(m)) {
#' look up neighbouring pixels
isna <- is.na(m)
rc <- nearest.valid.pixel(x[isna], y[isna], Zim, nsearch=2)
m[isna] <- Zim$v[cbind(rc$row, rc$col)]
}
}
)
return(m)
}
as.tess <- function(X) {
UseMethod("as.tess")
}
as.tess.tess <- function(X) {
fields <-
switch(X$type,
rect={ c("xgrid", "ygrid") },
tiled={ "tiles" },
image={ "image" },
stop(paste("Unrecognised tessellation type", sQuote(X$type))))
fields <- c(c("type", "window", "n", "marks"), fields)
X <- unclass(X)[fields]
class(X) <- c("tess", class(X))
return(X)
}
as.tess.im <- function(X) {
return(tess(image = X))
}
as.tess.list <- function(X) {
W <- lapply(X, as.owin)
return(tess(tiles=W))
}
as.tess.owin <- function(X) {
return(tess(tiles=list(X)))
}
domain.tess <- Window.tess <- function(X, ...) { as.owin(X) }
intersect.tess <- function(X, Y, ..., keepempty=FALSE, keepmarks=FALSE, sep="x") {
X <- as.tess(X)
check.1.string(sep)
if(is.owin(Y)) {
## intersection of a tessellation with a window
if(Y$type == "mask") {
## convert to pixel image
Xtiles <- tiles(X)
seqXtiles <- seq_along(Xtiles)
result <- as.im(Y, value=factor(1, levels=seqXtiles))
for(i in seqXtiles) {
tilei <- Xtiles[[i]]
result[tilei] <- i
}
result <- result[Y, drop=FALSE]
out <- tess(image=result, window=Y, keepempty=keepempty)
if(keepmarks && !is.null(marx <- marks(X))) {
if(keepempty) {
marks(out) <- marx
} else {
#' identify non-empty tiles
tab <- table(factor(result[], levels=seqXtiles))
marks(out) <- marksubset(marx, tab > 0)
}
}
return(out)
} else {
## efficient code when Y is a window, retaining names of tiles of X
Ztiles <- lapply(tiles(X), intersect.owin, B=Y, ..., fatal=FALSE)
isempty <- !keepempty & sapply(Ztiles, is.empty)
Ztiles <- Ztiles[!isempty]
Xwin <- as.owin(X)
Ywin <- Y
Zwin <- intersect.owin(Xwin, Ywin)
out <- tess(tiles=Ztiles, window=Zwin, keepempty=keepempty)
if(keepmarks) {
marx <- marks(X)
if(!is.null(marx))
marx <- as.data.frame(marx)[!isempty, , drop=FALSE]
marks(out) <- marx
}
return(out)
}
}
## general case: intersection of two tessellations
Y <- as.tess(Y)
Xtiles <- tiles(X)
Ytiles <- tiles(Y)
Ztiles <- list()
namesX <- tilenames(X)
namesY <- tilenames(Y)
if(keepmarks) {
## initialise the mark variables to be inherited from parent tessellations
Xmarks <- as.data.frame(marks(X))
Ymarks <- as.data.frame(marks(Y))
gotXmarks <- (ncol(Xmarks) > 0)
gotYmarks <- (ncol(Ymarks) > 0)
if(gotXmarks && gotYmarks) {
colnames(Xmarks) <- paste0("X", colnames(Xmarks))
colnames(Ymarks) <- paste0("Y", colnames(Ymarks))
}
if(gotXmarks || gotYmarks) {
marx <- if(gotXmarks && gotYmarks) {
cbind(Xmarks[integer(0), , drop=FALSE],
Ymarks[integer(0), , drop=FALSE])
} else if(gotXmarks) {
Xmarks[integer(0), , drop=FALSE]
} else {
Ymarks[integer(0), , drop=FALSE]
}
} else keepmarks <- FALSE
}
## now compute intersection tiles
Xtrivial <- (length(Xtiles) == 1)
for(i in seq_along(Xtiles)) {
Xi <- Xtiles[[i]]
Ti <- lapply(Ytiles, intersect.owin, B=Xi, ..., fatal=FALSE)
isempty <- !keepempty & sapply(Ti, is.empty)
nonempty <- !isempty
if(any(nonempty)) {
Ti <- Ti[nonempty]
names(Ti) <- if(Xtrivial) namesY[nonempty] else
paste(namesX[i], namesY[nonempty], sep=sep)
Ztiles <- append(Ztiles, Ti)
if(keepmarks) {
extra <- if(gotXmarks && gotYmarks) {
data.frame(X=Xmarks[i, ,drop=FALSE],
Y=Ymarks[nonempty, ,drop=FALSE],
row.names=NULL)
} else if(gotYmarks) {
Ymarks[nonempty, ,drop=FALSE]
} else {
Xmarks[rep(i, sum(nonempty)), ,drop=FALSE]
}
marx <- rbind(marx, extra)
}
}
}
## form tessellation object
Xwin <- as.owin(X)
Ywin <- as.owin(Y)
Zwin <- intersect.owin(Xwin, Ywin)
out <- tess(tiles=Ztiles, window=Zwin, keepempty=keepempty)
if(keepmarks)
marks(out) <- marx
return(out)
}
venn.tess <- function(..., window=NULL, labels=FALSE) {
argh <- list(...)
nargh <- length(argh)
if(nargh == 0) stop("No arguments given")
iswin <- sapply(argh, is.owin)
istes <- sapply(argh, is.tess)
if(!all(iswin | istes))
stop("All arguments must be windows or tessellations", call.=FALSE)
nama <- names(argh)
if(sum(nzchar(nama)) < nargh)
nama <- paste0("T", seq_len(nargh))
W <- window %orifnull% do.call(union.owin, unname(lapply(argh, as.owin)))
for(i in seq_len(nargh)) {
A <- argh[[i]]
if(is.owin(A)) {
Z <- list(A, Out=setminus.owin(W, A))
names(Z) <- paste0(c("", "not"), nama[i])
A <- tess(tiles=Z, window=W)
if(labels) marks(A) <- c(TRUE, FALSE)
}
Y <- if(i == 1) A else intersect.tess(Y, A, keepmarks=labels)
}
if(labels) colnames(marks(Y)) <- nama
return(Y)
}
bdist.tiles <- local({
vdist <- function(x,w) {
z <- as.ppp(vertices(x), W=w, check=FALSE)
min(bdist.points(z))
}
edist <- function(x,b) {
xd <- crossdist(edges(x, check=FALSE), b, type="separation")
min(xd)
}
bdist.tiles <- function(X) {
if(!is.tess(X))
stop("X must be a tessellation")
W <- as.owin(X)
switch(X$type,
rect=,
tiled={
tt <- tiles(X)
if(is.convex(W)) {
# distance is minimised at a tile vertex
d <- sapply(tt, vdist, w=W)
} else {
# coerce everything to polygons
W <- as.polygonal(W)
tt <- lapply(tt, as.polygonal)
# compute min dist from tile edges to window edges
d <- sapply(tt, edist, b=edges(W))
}
},
image={
Xim <- X$image
# compute boundary distance for each pixel
bd <- bdist.pixels(as.owin(Xim), style="image")
bd <- bd[W, drop=FALSE]
# split over tiles
bX <- split(bd, X)
# compute minimum distance over each level of factor
d <- sapply(bX, function(z) { summary(z)$min })
}
)
return(d)
}
bdist.tiles
})
## ......... geometrical transformations ..................
shift.tess <- function(X, ...) {
Y <- X
Y$window <- wY <- shift(X$window, ...)
vec <- getlastshift(wY)
switch(X$type,
rect={
Y$xgrid <- Y$xgrid + vec[1]
Y$ygrid <- Y$ygrid + vec[2]
},
tiled={
Y$tiles <- lapply(Y$tiles, shift, vec=vec)
},
image = {
Y$image <- shift(Y$image, vec)
})
attr(Y, "lastshift") <- vec
return(Y)
}
affine.tess <- function(X, mat=diag(c(1,1)), vec=c(0,0), ...) {
Y <- X
Y$window <- affine(X$window, mat=mat, vec=vec, ...)
switch(Y$type,
rect = {
if(all(mat == diag(diag(mat)))) {
## result is rectangular
Y$xgrid <- sort(mat[1,1] * X$xgrid + vec[1])
Y$ygrid <- sort(mat[2,2] * X$ygrid + vec[2])
} else {
## shear transformation; treat rectangles as general tiles
Y <- tess(tiles=tiles(X), window=Y$window)
Y$tiles <- lapply(Y$tiles, affine, mat=mat, vec=vec, ...)
}
},
tiled={
Y$tiles <- lapply(Y$tiles, affine, mat=mat, vec=vec, ...)
},
image = {
Y$image <- affine(Y$image, mat=mat, vec=vec, ...)
})
return(Y)
}
reflect.tess <- function(X) {
Y <- X
Y$window <- reflect(Y$window)
switch(X$type,
rect = {
Y$xgrid <- rev(- Y$xgrid)
Y$ygrid <- rev(- Y$ygrid)
},
tiled = {
Y$tiles <- lapply(Y$tiles, reflect)
},
image = {
Y$image <- reflect(Y$image)
})
return(Y)
}
flipxy.tess <- function(X) {
Y <- X
Y$window <- flipxy(Y$window)
switch(X$type,
rect = {
Y$xgrid <- X$ygrid
Y$ygrid <- X$xgrid
},
tiled = {
Y$tiles <- lapply(Y$tiles, flipxy)
},
image = {
Y$image <- flipxy(Y$image)
})
return(Y)
}
scalardilate.tess <- function(X, f, ...) {
Y <- X
Y$window <- scalardilate(X$window, f, ...)
switch(X$type,
rect = {
Y$xgrid <- f * Y$xgrid
Y$ygrid <- f * Y$ygrid
},
tiled = {
Y$tiles <- lapply(Y$tiles, scalardilate, f=f, ...)
},
image = {
Y$image <- scalardilate(Y$image, f=f, ...)
})
return(Y)
}
rotate.tess <- function(X, angle=pi/2, ..., centre=NULL) {
if(angle %% (2 * pi) == 0) return(X)
if(!is.null(centre)) {
X <- shift(X, origin=centre)
negorigin <- getlastshift(X)
} else negorigin <- NULL
Y <- X
Y$window <- rotate(X$window, angle=angle, ...)
switch(X$type,
rect = {
if(angle %% (pi/2) == 0) {
## result is rectangular
co <- round(cos(angle))
si <- round(sin(angle))
Y$xgrid <- sort((if(co == 0) 0 else (co * X$xgrid)) -
(if(si == 0) 0 else (si * X$ygrid)))
Y$ygrid <- sort((if(si == 0) 0 else (si * X$xgrid)) +
(if(co == 0) 0 else (co * X$ygrid)))
} else {
## general tessellation
Y <- tess(tiles=lapply(tiles(X), rotate, angle=angle, ...),
window=Y$window)
}
},
tiled = {
Y$tiles <- lapply(X$tiles, rotate, angle=angle, ...)
},
image = {
Y$image <- rotate(X$image, angle=angle, ...)
})
if(!is.null(negorigin))
Y <- shift(Y, -negorigin)
return(Y)
}
as.data.frame.tess <- function(x, ...) {
switch(x$type,
rect =,
tiled = {
y <- lapply(tiles(x), as.data.frame, ...)
z <- mapply(assignDFcolumn,
x=y, value=tilenames(x),
MoreArgs=list(name="Tile", ...),
SIMPLIFY=FALSE)
z <- do.call(rbind, z)
row.names(z) <- NULL
},
image = {
z <- as.data.frame(x$image, ...)
if(!is.na(m <- match("value", colnames(z))))
colnames(z)[m] <- "Tile"
},
{
z <- NULL
warning("Unrecognised type of tessellation")
})
return(z)
}
connected.tess <- function(X, ...) {
Xim <- as.im(X, ...)
X <- as.tess(Xim)
tilesX <- tiles(X)
namesX <- names(tilesX)
shards <- lapply(tilesX, connected) # list of factor images
shardnames <- lapply(shards, levels)
nshards <- lengths(shardnames)
broken <- (nshards > 1)
#' unbroken tiles keep their original tile names
shardnames[!broken] <- namesX[!broken]
#' shards of broken tiles are named "tilename[i] shard j"
shardnames[broken] <- mapply(paste,
namesX[broken], "shard", shardnames[broken],
SIMPLIFY=FALSE)
#' rename them
shards <- mapply("levels<-", shards, shardnames, SIMPLIFY=FALSE)
#' separate them
shards <- lapply(lapply(shards, as.tess), tiles)
shards <- unlist(shards, recursive=FALSE, use.names=FALSE)
names(shards) <- unlist(shardnames)
#' form tessellation
result <- tess(tiles=shards, window=as.owin(Xim))
result
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.