https://github.com/cran/spatstat
Raw File
Tip revision: 3aca716ce2576a0dab83f08052acd47afed8ee6a authored by Adrian Baddeley on 29 February 2012, 00:00:00 UTC
version 1.25-4
Tip revision: 3aca716
util.R
#
#    util.S    miscellaneous utilities
#
#    $Revision: 1.85 $    $Date: 2011/12/16 09:49:06 $
#
#  (a) for matrices only:
#
#    matrowany(X) is equivalent to apply(X, 1, any)
#    matrowall(X) "   "  " "  "  " apply(X, 1, all)
#    matcolany(X) "   "  " "  "  " apply(X, 2, any)
#    matcolall(X) "   "  " "  "  " apply(X, 2, all)
#
#  (b) for 3D arrays only:
#    apply23sum(X)  "  "   "  " apply(X, c(2,3), sum)
#
#  (c) weighted histogram
#    whist()
#
#  (d) for matrices:
#    matrixsample()
#           subsamples or supersamples a matrix
#
#
matrowsum <- function(x) {
  x %*% rep(1, ncol(x))
}

matcolsum <- function(x) {
  rep(1, nrow(x)) %*% x
}
  
matrowany <- function(x) {
  (matrowsum(x) > 0)
}

matrowall <- function(x) {
  (matrowsum(x) == ncol(x))
}

matcolany <- function(x) {
  (matcolsum(x) > 0)
}

matcolall <- function(x) {
  (matcolsum(x) == nrow(x))
}

########
    # hm, this is SLOWER

apply23sum <- function(x) {
  dimx <- dim(x)
  if(length(dimx) != 3)
    stop("x is not a 3D array")
  result <- array(0, dimx[-1])

  nz <- dimx[3]
  for(k in 1:nz) {
    result[,k] <- matcolsum(x[,,k])
  }
  result
}
    
#######################
#
#    whist      weighted histogram
#

whist <- function(x, breaks, weights=NULL) {
    N <- length(breaks)
    if(length(x) == 0) 
      h <- rep(0, N+1)
    else {
      # classify data into histogram cells (breaks need not span range of data)
      cell <- findInterval(x, breaks, rightmost.closed=TRUE)
      cell <- factor(cell, levels=0:N)
      # compute weighted histogram
      if(is.null(weights))
        h <- table(cell)
      else 
        h <- unlist(lapply(split(weights, cell), sum, na.rm=TRUE))
    }
    h <- as.numeric(h)
    y <- h[2:N]
    attr(y, "low") <- h[1]
    attr(y, "high") <- h[N+1]
    return(y)
}

######################
#
#   matrixsample         subsample or supersample a matrix
#

matrixsample <- function(mat, newdim, phase=c(0,0), scale, na.value=NA) {
  # 'phase+1' is the position of the [1,1] corner of the new matrix
  #  expressed in the coordinates of the old matrix.
  # 'scale' is the size of one step in the new matrix,
  #  expressed in the coordinates of the old matrix.
  # Both 'phase' and 'scale' can take any real value.
  olddim <- dim(mat)
  if(missing(scale)) scale <- (olddim - 1)/(newdim - 1)
  scale <- ensure2vector(scale)
  newdim  <- ensure2vector(newdim)
  newmat <- matrix(na.value, newdim[1], newdim[2])
  newrow <- 1:newdim[1]
  newcol <- 1:newdim[2]
  oldrow <- round(1 + phase[1] + (newrow-1) * scale[1])
  oldcol <- round(1 + phase[2] + (newcol-1) * scale[2])
  oldrow.ok <- (oldrow >= 1) & (oldrow <= olddim[1])
  oldcol.ok <- (oldcol >= 1) & (oldcol <= olddim[2])
  newmat[oldrow.ok, oldcol.ok] <- mat[oldrow[oldrow.ok],
                                      oldcol[oldcol.ok]]
  return(newmat)
}

# common invocation of matrixsample

rastersample <- function(X, Y) {
  stopifnot(is.im(X) || is.mask(X))
  stopifnot(is.im(Y) || is.mask(Y))
  phase <- c((Y$yrow[1] - X$yrow[1])/X$ystep,
             (Y$xcol[1] - X$xcol[1])/X$xstep)
  scale <- c(Y$ystep/X$ystep,
             Y$xstep/X$xstep)
  if(is.im(X)) {
    if(!is.im(Y)) Y <- as.im(Y)
    Y$v <- matrixsample(X$v, Y$dim, phase=phase, scale=scale, na.value=NA)
  } else {
    if(!is.mask(Y)) Y <- as.mask(Y)
    Y$m <- matrixsample(X$m, Y$dim, phase=phase, scale=scale, na.value=FALSE)
  }
  return(Y)
}

pointgrid <- function(W, ngrid) {
  W <- as.owin(W)
  masque <- as.mask(W, dimyx=ngrid)
  xx <- raster.x(masque)
  yy <- raster.y(masque)
  xx <- xx[masque$m]
  yy <- yy[masque$m]
  return(ppp(xx, yy, W))
}

# text magic

commasep <- function(x) {
  px <- paste(x)
  nx <- length(px)
  if(nx <= 1) return(px)
  commas <- c(rep(", ", length(px)-2),
              " and ", "")
  return(paste(paste(px, commas, sep=""), collapse=""))
}

paren <- function(x, type="(") {
  switch(type,
         "(" = {
           out <- paste("(", x, ")", sep="")
         },
         "[" = {
           out <- paste("[", x, "]", sep="")
         },
         "{" = {
           out <- paste("{", x, "}", sep="")
         })
  out
}

prange <- function(x) {
  stopifnot(length(x) == 2)
  paren(paste(x, collapse=", "), "[")
}

  
ordinal <- function(k) {
  last <- abs(k) %% 10
  lasttwo <- abs(k) %% 100
  isteen <- (lasttwo > 10 & lasttwo < 20)
  ending <- ifelse(isteen, "th",
                   ifelse(last == 1, "st",
                          ifelse(last == 2, "nd",
                                 ifelse(last == 3, "rd",
                                        "th"))))
  return(paste(k, ending, sep=""))
}

# equivalent to rev(cumsum(rev(x)))

revcumsum <- function(x) {
  n <- length(x)
  if(identical(storage.mode(x), "integer")) {
    z <- .C("irevcumsum",
            x=as.integer(x),
            as.integer(n),
            PACKAGE="spatstat")
    return(z$x)
  } else {
    z <- .C("drevcumsum",
            x=as.double(x),
            as.integer(n),
            PACKAGE="spatstat")
    return(z$x)
  }
}

prolongseq <- function(x, newrange) {
  stopifnot(length(newrange) == 2 && newrange[1] < newrange[2])
  stopifnot(length(x) >= 2)
  dx <- diff(x)
  if(any(dx <= 0))
    stop("x must be an increasing sequence")
  if(diff(range(dx)) > 0.01 * abs(mean(dx)))
    stop("x must be evenly spaced")
  dx <- mean(dx)

  # add or trim data to left
  if(x[1] > newrange[1]) {
    leftbit <- seq(from=x[1], to=newrange[1], by= -dx) 
    x <- c(rev(leftbit), x[-1])
  } else 
    x <- x[x >= newrange[1]]

  # add or trim data to right
  nx <- length(x)
  if(newrange[2] > x[nx]) {
    rightbit <- seq(from=x[nx], to=newrange[2], by= dx)
    x <- c(x[-nx], rightbit)
  } else 
    x <- x[x <= newrange[2]]

  return(x)
}

intersect.ranges <- function(a, b, fatal=TRUE) {
  lo <- max(a[1],b[1])
  hi <- min(a[2],b[2])
  if(lo >= hi) {
    if(fatal) stop("Intersection is empty")
    else return(NULL)
  }
  return(c(lo, hi))
}

assign(".Spatstat.ProgressBar", NULL, envir = .spEnv)

progressreport <- function(i, n, every=max(1, ceiling(n/100)),
                           nperline=min(charsperline,
                             every * ceiling(charsperline /(every+3))),
                           charsperline=60,
                           style=spatstat.options("progress")) {
  switch(style,
         txtbar={
           .Spatstat.ProgressBar <- get(".Spatstat.ProgressBar", envir = .spEnv)
           if(i == 1)
             assign(".Spatstat.ProgressBar", txtProgressBar(1, n, 1, style=3),
                    envir = .spEnv)
           setTxtProgressBar(.Spatstat.ProgressBar, i)
           if(i == n) {
             close(.Spatstat.ProgressBar)
             assign(".Spatstat.ProgressBar", NULL, envir = .spEnv)
           }
         },
         tty={
           if(i == n) 
             cat(paste(n, ".\n", sep=""))
           else if(every == 1 || i <= 3)
             cat(paste(i, ",", if(i %% nperline == 0) "\n" else " ", sep=""))
           else {
             if(i %% every == 0) 
               cat(i)
             else
               cat(".")
             if(i %% nperline == 0)
               cat("\n")
           }
           flush.console()
         },
         stop(paste("Unrecognised option for style:", dQuote(style)))
         )
  return(invisible(NULL))
}
  
numalign <- function(i, nmax, zero="0") {
  stopifnot(i <= nmax)
  nplaces <- as.integer(ceiling(log10(nmax+1)))
  out <- blank <- paste(rep(zero, nplaces), collapse="")
  istring <- paste(i)
  ilen <- nchar(istring)
  substr(out, nplaces-ilen+1, nplaces) <- istring
  return(out)
}

ensure2vector <- function(x) {
  xname <- deparse(substitute(x))
  if(!is.numeric(x))
    stop(paste(xname, "is not numeric"))
  n <- length(x)
  if(n == 0 || n > 2)
    stop(paste(xname, "should be of length 1 or 2"))
  if(n == 1)
    return(rep(x,2))
  return(x)
}

check.nvector <- function(v, npoints, fatal=TRUE, things="data points",
                          naok=FALSE) {
  # vector of numeric values for each point/thing
  vname <- sQuote(deparse(substitute(v)))
  whinge <- NULL
  if(!is.numeric(v))
    whinge <- paste(vname, "is not numeric")
  else if(length(v) != npoints)
    whinge <- paste("The length of", vname,
                    "should equal the number of", things)
  else if(!naok && any(is.na(v)))
    whinge <- paste("Some values of", vname, "are NA or NaN")
  #
  if(!is.null(whinge)) {
    if(fatal) stop(whinge) else return(FALSE)
  }
  return(TRUE)
}

check.nmatrix <- function(m, npoints, fatal=TRUE, things="data points",
                          naok=FALSE, squarematrix=TRUE) {
  # matrix of values for each thing or each pair of things
  mname <- sQuote(deparse(substitute(m)))
  whinge <- NULL
  if(!is.matrix(m))
    whinge <- paste(mname, "should be a matrix")
  else if(squarematrix && (nrow(m) != ncol(m)))
    whinge <- paste(mname, "should be a square matrix")
  else if(!naok && any(is.na(m)))
    whinge <- paste("Some values of", mname, "are NA or NaN")
  else if(nrow(m) != npoints)
    whinge <- paste("Number of rows in", mname,
               "does not match number of", things)
  #
  if(!is.null(whinge)) {
    if(fatal) stop(whinge) else return(FALSE)
  }
  return(TRUE)
}

check.named.vector <- function(x, nam, context="", namopt=character(0)) {
  xtitle <- deparse(substitute(x))
  check.named.thing(x, nam, namopt, sQuote(xtitle),
                    is.numeric(x), "vector", context)
  opt <- namopt %in% names(x)
  return(x[c(nam, namopt[opt])])
}

check.named.list <- function(x, nam, context="", namopt=character(0)) {
  xtitle <- deparse(substitute(x))
  check.named.thing(x, nam, namopt, sQuote(xtitle),
                    is.list(x), "list", context)  
  opt <- namopt %in% names(x)
  return(x[c(nam, namopt[opt])])
}

check.named.thing <- function(x, nam, namopt=character(0), xtitle=NULL,
                              valid=TRUE, type="object", context="",
                              fatal=TRUE) {
  if(is.null(xtitle))
    xtitle <- sQuote(deparse(substitute(x)))
  # check whether names(x) contains all obligatory names 'nam'
  # and possibly some of the optional names 'namopt'
  namesx <- names(x)
  omitted <- !(nam %in% namesx)
  foreign <- !(namesx %in% c(nam, namopt))
  if(valid && !any(omitted) && !any(foreign))
    return(character(0))
  # some condition violated
  if(nzchar(context))
    xtitle <- paste(context, xtitle)
  whinge <- paste(xtitle,
                  "must be a named", paste(type, ",", sep=""),
                  "with components", commasep(nam))
  if(length(namopt) > 0)
    whinge <- paste(whinge, paren(paste("and optionally", commasep(namopt))))
  if(any(omitted)) {
    grizzle <- paste(ngettext(sum(omitted), "parameter", "parameters"),
                     commasep(nam[omitted]),
                     "omitted")
    whinge <- paste(whinge, grizzle, sep="; ")
  }
  if(any(foreign)) {
    grizzle <- paste(ngettext(sum(foreign), "component", "components"),
                     commasep(namesx[foreign]),
                     "not recognised")
    whinge <- paste(whinge, grizzle, sep="; ")
  }
  if(fatal)
    stop(whinge, call.=FALSE)
  return(whinge)
}

forbidNA <- function(x, context="", xname) {
  if(missing(xname)) xname <- deparse(substitute(x))
  if(any(is.na(x))) {
    offence <- ngettext(length(x), "be NA", "contain NA values")
    whinge <- paste(context, sQuote(xname), "must not", offence)
    stop(whinge, call.=FALSE)
  }
  return(TRUE)
}

check.finite <- function(x, context="", xname) {
  if(missing(xname)) xname <- deparse(substitute(x))
  forbidNA(x, context, xname)
  if(any(!is.finite(x))) {
    oblige <- ngettext(length(x), "be a finite value", "contain finite values")
    whinge <- paste(context, sQuote(xname), "must", oblige)
    stop(whinge, call.=FALSE)
  }
  return(TRUE)
}

evenly.spaced <- function(x, tol=1e-07) {
  # test whether x is evenly spaced and increasing
  dx <- diff(x)
  if(any(dx <= .Machine$double.eps))
    return(FALSE)
  # The following test for equal spacing is used in hist.default
  if(diff(range(dx)) > tol * mean(dx))
    return(FALSE)
  return(TRUE)
}

adjustthinrange <- function(ur,vstep,vr) {
  if(diff(ur) >= vstep) return(ur)
  ur <- mean(ur) + c(-1,1) * vstep/2
  if(ur[1] < vr[1]) ur <- vr[1] + c(0,1)*vstep
  if(ur[2] > vr[2]) ur <- vr[2] - c(1,0)*vstep
  return(ur)
}

validposint <- function(n, caller, fatal=TRUE) {
  nname <- deparse(substitute(n))
  if(length(n) != 1 || n != round(n) || n <=0) {
    if(!fatal)
      return(FALSE)
    prefix <- if(!missing(caller)) paste("In ", caller, ",", sep="") else NULL
    stop(paste(prefix, nname, "should be a single positive integer"),
         call.=FALSE)
  }
  return(TRUE)
}

# wrangle data.frames

firstfactor <- function(x) {
  stopifnot(is.data.frame(x))
  isfac <- unlist(lapply(as.list(x), is.factor))
  if(!any(isfac)) 
    return(NULL)
  return(x[, min(which(isfac))])
}

onecolumn <- function(m) {
  switch(markformat(m),
         none=stop("No marks provided"),
         vector=m,
         dataframe=m[,1, drop=TRUE],
         NA)
}

# errors and checks

complaining <- function(whinge, fatal=FALSE, value=NULL) {
  if(fatal) stop(whinge, call.=FALSE)
  warning(whinge, call.=FALSE)
  return(value)
}

check.1.real <- function(x, context="", fatal=TRUE) {
  xname <- deparse(substitute(x))
  if(!is.numeric(x) || length(x) != 1) {
    whinge <-  paste(xname, "should be a single number")
    if(nzchar(context)) whinge <- paste(context, whinge)
    return(complaining(whinge, fatal=fatal, value=FALSE))
  }
  return(TRUE)
}
   
explain.ifnot <- function(expr, context="") {
  ex <- deparse(substitute(expr))
  ans <- expr
  if(!identical(ans, TRUE))
    stop(paste(context, "it must be TRUE that", sQuote(ex)), call.=FALSE)
}

warn.ignored.args <- function(..., context=NULL) {
  if((narg <- length(list(...))) > 0) {
    whinge <- paste(narg, "unrecognised",
                    ngettext(narg, "argument was", "arguments were"),
                    "ignored")
    if(!is.null(context)) whinge <- paste(context, whinge)
    warning(context)
  }
}

multiply.only.finite.entries <- function(x, a) {
  # In ppm a potential value that is -Inf must remain -Inf
  # and a potential value that is 0 multiplied by NA remains 0
  ifelse(is.finite(x) & x != 0, x * a, x)
}

singlestring <- function(s, coll="") {
  s <- as.character(s)
  if(length(s) > 1)
    s <- paste(s, collapse=coll)
  return(s)
}

verbalogic <- function(x, op="and") {
  stopifnot(is.character(x))
  istrue <- (x == "TRUE")
  isfalse <- (x == "FALSE")
  isvariable <- !istrue & !isfalse
  y <- x[isvariable]
  switch(op,
         and={
           if(any(isfalse))
             return("FALSE")
           if(all(istrue))
             return("TRUE")
           return(paste(y, collapse=" and "))
         },
         or={
           if(all(isfalse))
             return("FALSE")
           if(any(istrue))
             return("TRUE")
           return(paste(y, collapse=" or "))
         },
         not={
           x[isfalse] <- "TRUE"
           x[istrue] <- "FALSE"
           x[isvariable] <- paste("not {", y, "}")
         },
         stop(paste("Unrecognised operation", sQuote(op))))
}

sensiblevarname <- function(guess, fallback, maxlen=12) {
  out <- if(is.character(guess) &&
            length(guess) == 1  &&
            make.names(guess) == guess) guess else fallback
  out <- substr(out, 1, maxlen)
  return(out)
}

good.names <- function(nama, defaults, suffices) {
  # ensure sensible, unique names 
  stopifnot(is.character(defaults))
  if(!missing(suffices))
    defaults <- paste(defaults, suffices, sep="")
  result <- nama
  if(is.null(result))
    result <- defaults
  else if(any(blank <- !nzchar(result)))
    result[blank] <- defaults[blank]
  if(any(duplicated(result)))
    result <- make.names(result, unique=TRUE)
  return(result)
}

    
cat.factor <- function (..., recursive=FALSE) {
  lll <- list(...)
  chk <- sapply(lll,is.factor)
  if(!all(chk))
    stop("First argument is a factor and at least one other argument is not.\n")
  lll <- lapply(lll,as.data.frame,nm="v1")
  return(do.call(rbind,lll)[,1])
}

nzpaste <- function(..., sep=" ", collapse=NULL) {
  # Paste only the non-empty strings
  v <- list(...)
  ok <- unlist(lapply(v, function(z) {any(nzchar(z))}))
  do.call("paste", append(v[ok], list(sep=sep, collapse=collapse)))
}

is.parseable <- function(x) {
  unlist(lapply(x, function(z) {
    !inherits(try(parse(text=z), silent=TRUE), "try-error")
  }))
}

make.parseable <- function(x) {
  if(all(is.parseable(x))) x else make.names(x)
}

# paste(expression(..)) seems to be broken

paste.expr <- function(x) {
  unlist(lapply(x, function(z) { paste(deparse(z), collapse="") }))
}

badprobability <- function(x, NAvalue=NA) {
  ifelse(is.na(x), NAvalue, !is.finite(x) | x < 0 | x > 1)
}

# test for equivalence of two functions 
samefunction <- function(f, g) {
  identical(deparse(f), deparse(g))
}

codetime <- local({
  uname <- c("min", "hours", "days", "weeks", "years",
             "thousand years", "million years", "billion years")
  multiple <- c(60, 60, 24, 7, 365/7, 1e3, 1e3, 1e3)
  function(x) {
    u <- "sec"
    for(k in seq_along(multiple)) {
      if(x > multiple[k]) {
        x <- x/multiple[k]
        u <- uname[k]
      } else break
    }
    paste(round(x, 1), u)
  }
})
back to top