# # ppp.S # # A class 'ppp' to define point patterns # observed in arbitrary windows in two dimensions. # # $Revision: 4.73 $ $Date: 2010/10/28 02:08:37 $ # # A point pattern contains the following entries: # # $window: an object of class 'owin' # defining the observation window # # $n: the number of points (for efficiency) # # $x: # $y: vectors of length n giving the Cartesian # coordinates of the points. # # It may also contain the entry: # # $marks: a vector of length n # whose entries are interpreted as the # 'marks' attached to the corresponding points. # #-------------------------------------------------------------------------- ppp <- function(x, y, ..., window, marks, check=TRUE ) { # Constructs an object of class 'ppp' # if(!missing(window)) verifyclass(window, "owin") else window <- owin(...) if(missing(x) && missing(y)) x <- y <- numeric(0) n <- length(x) if(length(y) != n) stop("coordinate vectors x and y are not of equal length") # validate x, y coordinates stopifnot(is.numeric(x)) stopifnot(is.numeric(y)) ok <- is.finite(x) & is.finite(y) if(any(!ok)) { nbg <- is.na(x) | is.na(y) if(any(nbg)) { howmany <- if(all(nbg)) "all" else paste(sum(nbg), "out of", length(nbg)) stop(paste(howmany, "coordinate values are NA or NaN")) } howmany <- if(!any(ok)) "all" else paste(sum(!ok), "out of", length(ok)) stop(paste(howmany, "coordinate values are infinite")) } names(x) <- NULL names(y) <- NULL # check (x,y) points lie inside window if(check && n > 0) { ok <- inside.owin(x, y, window) nout <- sum(!ok) if(nout > 0) { warning(paste(nout, ngettext(nout, "point was", "points were"), "rejected as lying outside the specified window")) rr <- ripras(x,y) bb <- bounding.box.xy(x,y) bb <- bounding.box(rr, bb, window) rejectwindow <- if(!is.null(rr)) rebound.owin(rr, bb) else bb rejects <- ppp(x[!ok], y[!ok], window=rejectwindow, check=FALSE) # discard illegal points x <- x[ok] y <- y[ok] n <- length(x) } } else nout <- 0 # initialise ppp object pp <- list(window=window, n=n, x=x, y=y) # coerce marks to appropriate forma if(missing(marks)) marks <- NULL if(is.hyperframe(marks)) stop("Hyperframes of marks are not implemented for ppp objects; use ppx") if(is.matrix(marks)) marks <- as.data.frame(marks) if(is.data.frame(marks)) { nc <- ncol(marks) if(nc == 0) marks <- NULL else if(nc == 1) marks <- marks[,,drop=TRUE] } # attach marks if(is.null(marks)) { # no marks pp$markformat <- "none" } else if(is.data.frame(marks)) { # data frame of marks pp$markformat <- "dataframe" if(nout > 0) { marks <- marks[ok, ] marks(rejects) <- marks[!ok,] } if(nrow(marks) != n) stop("number of rows of marks != length of x and y") pp$marks <- marks } else { # should be a vector or factor # To recognise vector, strip attributes if(!is.factor(marks)) attributes(marks) <- NULL if(!(is.vector(marks) || is.factor(marks))) stop("Format of marks not understood") # OK, it's a vector or factor pp$markformat <- "vector" if(nout > 0) { marks(rejects) <- marks[!ok] marks <- marks[ok] } if(length(marks) != n) stop("length of marks vector != length of x and y") names(marks) <- NULL pp$marks <- marks } class(pp) <- "ppp" if(check && any(duplicated(pp))) warning("data contain duplicated points") if(nout > 0) attr(pp, "rejects") <- rejects pp } # #-------------------------------------------------------------------------- # is.ppp <- function(x) { inherits(x, "ppp") } # #-------------------------------------------------------------------------- # as.ppp <- function(X, ..., fatal=TRUE) { UseMethod("as.ppp") } as.ppp.ppp <- function(X, ..., fatal=TRUE) { check <- resolve.defaults(list(...), list(check=FALSE))$check return(ppp(X$x, X$y, window=X$window, marks=X$marks, check=check)) } as.ppp.quad <- function(X, ..., fatal=TRUE) { return(union.quad(X)) } as.ppp.data.frame <- function(X, W = NULL, ..., fatal=TRUE) { check <- resolve.defaults(list(...), list(check=TRUE))$check if(!(ncol(X) %in% 2:3)) { if(fatal) stop("X must be a two-column or three-column data frame") else return(NULL) } if(is.null(W)) { if(fatal) stop("x,y coords given but no window specified") else return(NULL) } if(is.function(W)) Z <- cobble.xy(X[,1], X[,2], W, fatal) else { win <- as.owin(W) Z <- ppp(X[,1], X[,2], window = win, check=check) } # add marks from third column if(ncol(X) == 3) Z <- Z %mark% X[,3] return(Z) } as.ppp.matrix <- function(X, W = NULL, ..., fatal=TRUE) { check <- resolve.defaults(list(...), list(check=TRUE))$check if(!verifyclass(X, "matrix", fatal=fatal) || !is.numeric(X)) { if(fatal) stop("X must be a two-column matrix of numbers") else return(NULL) } if(!(ncol(X) %in% 2:3)) { if(fatal) stop("X must be a two-column or three-column matrix of numbers") else return(NULL) } if(is.null(W)) { if(fatal) stop("x,y coords given but no window specified") else return(NULL) } if(is.function(W)) Z <- cobble.xy(X[,1], X[,2], W, fatal) else { win <- as.owin(W) Z <- ppp(X[,1], X[,2], window = win, check=check) } # add marks from third column if(ncol(X) == 3) Z <- Z %mark% X[,3] return(Z) } as.ppp.default <- function(X, W=NULL, ..., fatal=TRUE) { # tries to coerce data X to a point pattern # X may be: # 1. a structure with entries x, y, xl, xu, yl, yu # 2. a structure with entries x, y, area where # 'area' has entries xl, xu, yl, yu # 3. a structure with entries x, y # 4. a vector of length 2, interpreted as a single point. # The second argument W is coerced to an object of class 'owin' by the # function "as.owin" in window.S # If X also has an entry X$marks # then this will be interpreted as the marks vector for the pattern. # check <- resolve.defaults(list(...), list(check=TRUE))$check if(checkfields(X, c("x", "y", "xl", "xu", "yl", "yu"))) { xrange <- c(X$xl, X$xu) yrange <- c(X$yl, X$yu) if(is.null(X$marks)) Z <- ppp(X$x, X$y, xrange, yrange, check=check) else Z <- ppp(X$x, X$y, xrange, yrange, marks=X$marks, check=check) return(Z) } else if(checkfields(X, c("x", "y", "area")) && checkfields(X$area, c("xl", "xu", "yl", "yu"))) { win <- as.owin(X$area) if (is.null(X$marks)) Z <- ppp(X$x, X$y, window=win, check=check) else Z <- ppp(X$x, X$y, window=win, marks = X$marks, check=check) return(Z) } else if(checkfields(X, c("x", "y"))) { if(is.function(W)) return(cobble.xy(X$x, X$y, W, fatal)) if(is.null(W)) { if(fatal) stop("x,y coords given but no window specified") else return(NULL) } win <- as.owin(W) if(is.null(X$marks)) Z <- ppp(X$x, X$y, window=win, check=check) else Z <- ppp(X$x, X$y, window=win, marks=X$marks, check=check) return(Z) } else if(is.vector(X) && length(X) == 2) { win <- as.owin(W) Z <- ppp(X[1], X[2], window=win, check=check) return(Z) } else { if(fatal) stop("Can't interpret X as a point pattern") else return(NULL) } } cobble.xy <- function(x, y, f=ripras, fatal=TRUE) { if(!is.function(f)) stop("f is not a function") w <- f(x,y) if(!is.owin(w)) { gripe <- "Supplied function f did not return an owin object" if(fatal) stop(gripe) else { warning(gripe) return(NULL) } } return(ppp(x, y, window=w)) } # -------------------------------------------------------------- "[.ppp" <- function(x, i, j, drop, ...) { verifyclass(x, "ppp") if(missing(i) && missing(j)) return(x) if(!missing(i)) { if(inherits(i, "owin")) { # i is a window window <- i ok <- inside.owin(x$x, x$y, window) x <- ppp(x$x[ok], x$y[ok], window=window, #SIC marks=marksubset(x$marks, ok), check=FALSE) } else if(inherits(i, "im")) { # i is an image if(i$type != "logical") stop(paste("Subset operator X[i] undefined", "when i is a pixel image", "unless it has logical values"), call.=FALSE) # convert logical image to window e <- sys.frame(sys.nframe()) window <- solutionset(i, e) ok <- inside.owin(x$x, x$y, window) x <- ppp(x$x[ok], x$y[ok], window=window, #SIC marks=marksubset(x$marks, ok), check=FALSE) } else { # assume i is a subset index if(x$n == 0) return(x) subset <- i x <- ppp(x$x[subset], x$y[subset], window=x$window, marks=marksubset(x$marks, subset), check=FALSE) } } if(!missing(j)) x <- x[j] # invokes code above return(x) } # ------------------------------------------------------------------ # # scanpp <- function(filename, window, header=TRUE, dir="", multitype=FALSE) { filename <- if(dir=="") filename else paste(dir, filename, sep=.Platform$file.sep) df <- read.table(filename, header=header) if(header) { # check whether there are columns named 'x' and 'y' colnames <- dimnames(df)[[2]] xycolumns <- match(c("x", "y"), colnames, 0) named <- all(xycolumns > 0) } else { named <- FALSE } if(named) { x <- df$x y <- df$y } else { # assume x, y given in columns 1, 2 respectively x <- df[,1] y <- df[,2] xycolumns <- c(1,2) } if(ncol(df) == 2) X <- ppp(x, y, window=window) else { marks <- df[ , -xycolumns] if(multitype) marks <- factor(marks) X <- ppp(x, y, window=window, marks = marks) } X } #------------------------------------------------------------------- "markspace.integral" <- function(X) { verifyclass(X, "ppp") if(!is.marked(X, dfok=TRUE)) return(1) if(is.multitype(X)) return(length(levels(marks(X)))) else stop("Don't know how to compute total mass of mark space") } #------------------------------------------------------------------- print.ppp <- function(x, ...) { verifyclass(x, "ppp") ism <- is.marked(x, dfok=TRUE) cat(paste(if(ism) "marked" else NULL, "planar point pattern:", x$n, ngettext(x$n, "point", "points"), "\n")) if(ism) { mks <- marks(x, dfok=TRUE) if(is.data.frame(mks)) { # data frame of marks cat(paste("Mark variables: ", paste(names(mks), collapse=", "), "\n")) } else { # vector of marks if(is.factor(mks)) { cat("multitype, with ") cat(paste("levels =", paste(levels(mks), collapse="\t"),"\n")) } else { cat(paste("marks are", if(is.numeric(mks)) "numeric,", "of type", sQuote(typeof(mks)), "\n")) } } } print(x$window) if(!is.null(rejects <- attr(x, "rejects"))) { nrejects <- rejects$n cat(paste("\n***", nrejects, ngettext(nrejects, "illegal point", "illegal points"), "stored in", paste("attr(,", dQuote("rejects"), ")", sep=""), "***\n")) } return(invisible(NULL)) } summary.ppp <- function(object, ..., checkdup=TRUE) { verifyclass(object, "ppp") result <- list() result$is.marked <- is.marked(object, dfok=TRUE) result$n <- object$n result$window <- summary(object$window) result$intensity <- result$n/result$window$area if(checkdup) result$nduplicated <- sum(duplicated(object)) if(result$is.marked) { mks <- marks(object, dfok=TRUE) if(result$multiple.marks <- is.data.frame(mks)) { result$marknames <- names(mks) result$is.numeric <- FALSE result$marktype <- "dataframe" result$is.multitype <- FALSE } else { result$is.numeric <- is.numeric(mks) result$marknames <- "marks" result$marktype <- typeof(mks) result$is.multitype <- is.multitype(object) } if(result$is.multitype) { tm <- as.vector(table(mks)) tfp <- data.frame(frequency=tm, proportion=tm/sum(tm), intensity=tm/result$window$area, row.names=levels(mks)) result$marks <- tfp } else result$marks <- summary(mks) } class(result) <- "summary.ppp" if(!is.null(rejects <- attr(object, "rejects"))) result$rejects <- rejects$n return(result) } print.summary.ppp <- function(x, ..., dp=3) { verifyclass(x, "summary.ppp") cat(paste(if(x$is.marked) "Marked planar " else "Planar ", "point pattern: ", x$n, " points\n", sep="")) oneline <- resolve.defaults(list(...), list(oneline=FALSE))$oneline if(oneline) return(invisible(NULL)) unitinfo <- summary(x$window$units) cat(paste("Average intensity", signif(x$intensity,dp), "points per square", unitinfo$singular, unitinfo$explain, "\n")) ndup <- x$nduplicated if((!is.null(ndup)) && (ndup > 0)) cat("\n*Pattern contains duplicated points*\n") if(x$is.marked) { if(x$multiple.marks) { cat(paste("Mark variables: ", paste(x$marknames, collapse=", "), "\n")) cat("Summary:\n") print(x$marks) } else if(x$is.multitype) { cat("Multitype:\n") print(signif(x$marks,dp)) } else { cat(paste("marks are ", if(x$is.numeric) "numeric, ", "of type", sQuote(x$marktype), "\n")) cat("Summary:\n") print(x$marks) } } cat("\n") print(x$window) if(!is.null(nrejects <- x$rejects)) cat(paste("\n***", nrejects, ngettext(nrejects, "illegal point", "illegal points"), "stored in", paste("attr(,", dQuote("rejects"), ")", sep=""), "***\n")) return(invisible(x)) } # --------------------------------------------------------------- identify.ppp <- function(x, ...) { verifyclass(x, "ppp") if(!is.marked(x) || "labels" %in% names(list(...))) identify(x$x, x$y, ...) else { marx <- marks(x, dfok=FALSE) marques <- if(is.numeric(marx)) paste(signif(marx, 3)) else paste(marx) id <- identify(x$x, x$y, labels=marques, ...) mk <- marx[id] if(is.factor(marx)) mk <- levels(marx)[mk] cbind(id=id, marks=mk) } } rebound <- function(x, rect) { UseMethod("rebound") } rebound.ppp <- function(x, rect) { verifyclass(x, "ppp") x$window <- rebound.owin(x$window, rect) return(x) } as.data.frame.ppp <- function(x, row.names=NULL, ...) { df <- data.frame(x=x$x, y=x$y, row.names=row.names) marx <- marks(x, dfok=TRUE) if(is.null(marx)) return(df) if(is.data.frame(marx)) df <- cbind(df, marx) else df <- data.frame(df, marks=marx) return(df) } is.empty.ppp <- function(x) { return(x$n == 0) } npoints <- function(x) { UseMethod("npoints") } npoints.ppp <- function(x) { x$n }