R Under development (unstable) (2023-07-19 r84711 ucrt) -- "Unsuffered Consequences" Copyright (C) 2023 The R Foundation for Statistical Computing Platform: x86_64-w64-mingw32/x64 R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > #' > #' Header for all (concatenated) test files > #' > #' Require spatstat.geom > #' Obtain environment variable controlling tests. > #' > #' $Revision: 1.5 $ $Date: 2020/04/30 05:31:37 $ > > require(spatstat.geom) Loading required package: spatstat.geom Loading required package: spatstat.data spatstat.geom 3.2-4 > FULLTEST <- (nchar(Sys.getenv("SPATSTAT_TEST", unset="")) > 0) > ALWAYS <- TRUE > cat(paste("--------- Executing", + if(FULLTEST) "** ALL **" else "**RESTRICTED** subset of", + "test code -----------\n")) --------- Executing **RESTRICTED** subset of test code ----------- > # > # tests/nndist.R > # > # Check that nndist and nnwhich give > # results consistent with direct calculation from pairdist > # > # Similarly for nncross and distfun > # > # Also test whether minnndist(X) == min(nndist(X)) > # > # $Revision: 1.39 $ $Date: 2021/05/20 09:31:23 $ > # > > > local({ + eps <- sqrt(.Machine$double.eps) + f <- function(mat,k) { apply(mat, 1, function(z,n) { sort(z)[n] }, n=k+1) } + g <- function(mat,k) { apply(mat, 1, function(z,n) { order(z)[n] }, n=k+1) } + + ## ....... Two dimensions ................ + if(ALWAYS) { + X <- runifrect(24) + + nn <- nndist(X) + nnP <- f(pairdist(X), 1) + if(any(abs(nn - nnP) > eps)) + stop("nndist.ppp does not agree with pairdist") + + nn5 <- nndist(X, k=5) + nn5P <- f(pairdist(X), 5) + if(any(abs(nn5 - nn5P) > eps)) + stop("nndist.ppp(k=5) does not agree with pairdist") + + nw <- nnwhich(X) + nwP <- g(pairdist(X), 1) + if(any(nw != nwP)) + stop("nnwhich.ppp does not agree with pairdist") + + nw5 <- nnwhich(X, k=5) + nw5P <- g(pairdist(X), 5) + if(any(nw5 != nw5P)) + stop("nnwhich.ppp(k=5) does not agree with pairdist") + } + + if(FULLTEST) { + a <- nndist(X, method="test") + b <- nnwhich(X, method="test") + a <- nndist(X, method="test", k=1:2) + b <- nnwhich(X, method="test", k=1:2) + a2 <- nndist(cells[1:3], k=1:3) + b2 <- nnwhich(cells[1:3], k=1:3) + a3 <- nndist(cells[1]) + b3 <- nnwhich(cells[1]) + m <- factor((1:npoints(X)) %% 2 == 0) + a4 <- nndist.default(X, by=m, k=2) + b4 <- nnwhich.default(X, by=m, k=2) + } + + if(ALWAYS) { + ## nncross.ppp without options + Y <- runifrect(30) + Y <- Y[nndist(Y) > 0.02] + nc <- nncross(X,Y) + ncd <- nc$dist + ncw <- nc$which + cd <- crossdist(X,Y) + cdd <- apply(cd, 1, min) + cdw <- apply(cd, 1, which.min) + if(any(abs(ncd - cdd) > eps)) + stop("nncross()$dist does not agree with apply(crossdist(), 1, min)") + if(any(ncw != cdw)) + stop("nncross()$which does not agree with apply(crossdist(), 1, which.min)") + + ## nncross with sort on x + nc <- nncross(X,Y, sortby="x") + ncd <- nc$dist + ncw <- nc$which + if(any(abs(ncd - cdd) > eps)) + stop("nncross(sortby=x)$dist does not agree with apply(crossdist(), 1, min)") + if(any(ncw != cdw)) + stop("nncross(sortby=x)$which does not agree with apply(crossdist(), 1, which.min)") + + ## nncross with data pre-sorted on x + Y <- Y[order(Y$x)] + nc <- nncross(X,Y, is.sorted.Y=TRUE, sortby="x") + ncd <- nc$dist + ncw <- nc$which + cd <- crossdist(X,Y) + cdd <- apply(cd, 1, min) + cdw <- apply(cd, 1, which.min) + if(any(abs(ncd - cdd) > eps)) + stop("For sorted data, nncross()$dist does not agree with apply(crossdist(), 1, min)") + if(any(ncw != cdw)) + stop("For sorted data, nncross()$which does not agree with apply(crossdist(), 1, which.min)") + + ## sanity check for nncross with k > 1 + ndw <- nncross(X, Y, k=1:4, what="which") + if(any(is.na(ndw))) + stop("NA's returned by nncross.ppp(k > 1, what='which')") + nnc4 <- nncross(X, Y, k=1:4) + iswhich <- (substr(colnames(nnc4), 1, nchar("which")) == "which") + ndw <- nnc4[,iswhich] + if(any(is.na(ndw))) + stop("NA's returned by nncross.ppp(k > 1)$which") + + ## test of correctness for nncross with k > 1 + flipcells <- flipxy(cells) + calcwhich <- nncross(cells, flipcells, k=1:4, what="which") + truewhich <- t(apply(crossdist(cells,flipcells), 1, order))[,1:4] + if(any(calcwhich != truewhich)) + stop("nncross(k > 1) gives wrong answer") + } + + if(FULLTEST) { + ## example from Hank Stevens + A <- data.frame( + m= c("K", "K", "A1", "A2", "G", "A2", "A3"), + x=c(4.85, 6.76, 10.58, 19.18, 15.74, 19.08, 12.27), + y=c(5.60, 12.92, 11.14, 17.22, 5.74, 1.24, 2.20), + stringsAsFactors=TRUE + ) + X <- with(A, ppp(x, y, marks=m, window=bounding.box.xy(x, y))) + suspect <- nncross(X, X[7], iX=1:7, iY=7L)$dist + correct <- c(pairdist(X)[1:6, 7], Inf) + maxer <- max(abs(suspect[1:6] - correct[1:6])) + if(maxer > 0.001) + stop("Error in nncross (Inf values) in Hank Stevens example") + if(suspect[7] != Inf) + stop("Error in nncross (finite values) in Hank Stevens example") + M <- as.matrix(minnndist(X, by=marks(X))) + M[is.infinite(M)] <- 0 + maxer <- max(abs(M - t(M))) + if(maxer > 0.001) + stop("Error in minnndist(by) in Hank Stevens example") + } + + if(ALWAYS) { + #' cover some C code blocks + Z <- runifrect(50) + X <- Z[1:30] + Y <- Z[20:50] + iX <- 1:30 + iY <- 20:50 + Ndw <- nncross(X,Y, iX, iY, k=3) + Nw <- nncross(X,Y, iX, iY, k=3, what="which") + Nd <- nncross(X,Y, iX, iY, k=3, what="dist") + } + + if(FULLTEST) { + ## special cases + nndist(X[FALSE]) + nndist(X[1]) + nndist(X[1:3], k=4) + nndist(X[1:3], k=1:4) + nnwhich(X[FALSE]) + nnwhich(X[1]) + nnwhich(X[1:3], k=4) + nnwhich(X[1:3], k=1:4) + nncross(X[1:3], Y[FALSE]) + nncross(X[1:3], Y[1]) + nncross(X[1:3], Y[1:3], k=4) + nncross(X[1:3], Y[1:3], k=1:4) + } + + ## ....... Three dimensions ................ + + if(ALWAYS) { + rthree <- function(n) { pp3(runif(n), runif(n), runif(n), box3(c(0,1))) } + XX <- rthree(42) + X <- XX[1:20] + nn <- nndist(X) + nnP <- f(pairdist(X), 1) + if(any(abs(nn - nnP) > eps)) + stop("nndist.pp3 does not agree with pairdist") + + nn5 <- nndist(X, k=5) + nn5P <- f(pairdist(X), 5) + if(any(abs(nn5 - nn5P) > eps)) + stop("nndist.pp3(k=5) does not agree with pairdist") + + nw <- nnwhich(X) + nwP <- g(pairdist(X), 1) + if(any(nw != nwP)) + stop("nnwhich.pp3 does not agree with pairdist") + + nw5 <- nnwhich(X, k=5) + nw5P <- g(pairdist(X), 5) + if(any(nw5 != nw5P)) + stop("nnwhich.pp3(k=5) does not agree with pairdist") + + ff <- function(mat,k) { apply(mat, 1, function(z,n) { sort(z)[n] }, n=k) } + gg <- function(mat,k) { apply(mat, 1, function(z,n) { order(z)[n] }, n=k) } + + Y <- rthree(20) + Y <- Y[nndist(Y) > 0.02] + DXY <- crossdist(X,Y) + a <- nncross(X,Y) + a <- nncross(X,Y, what="dist") + a <- nncross(X,Y, what="which") + if(any(a != gg(DXY, 1))) + stop("incorrect result from nncross.pp3(what='which')") + a2 <- nncross(X,Y, k=2) + a2 <- nncross(X,Y, what="dist", k=2) + a2 <- nncross(X,Y, what="which", k=2) + if(any(a2 != gg(DXY, 2))) + stop("incorrect result from nncross.pp3(k=2, what='which')") + } + + if(FULLTEST) { + X <- XX + iX <- 1:42 + iZ <- 30:42 + Z <- X[iZ] + b <- nncross(X, Z, iX=iX, iY=iZ) + b <- nncross(X, Z, iX=iX, iY=iZ, what="which") + b <- nncross(X, Z, iX=iX, iY=iZ, what="dist") + b2 <- nncross(X, Z, iX=iX, iY=iZ, k=2) + b2 <- nncross(X, Z, iX=iX, iY=iZ, what="which", k=2) + b2 <- nncross(X, Z, iX=iX, iY=iZ, what="dist", k=2) + e1 <- nncross(X, Y[1:3], k=2:4) + c1 <- nncross(X, Y, sortby="var") + c2 <- nncross(X, Y, sortby="x") + c3 <- nncross(X, Y, sortby="y") + c4 <- nncross(X, Y, sortby="z") + Xsort <- X[order(coords(X)$x)] + c5 <- nncross(Xsort, Y, is.sorted.X=TRUE, sortby="x") + Ysort <- Y[order(coords(Y)$x)] + c6 <- nncross(Xsort, Ysort, is.sorted.X=TRUE, is.sorted.Y=TRUE, sortby="x") + } + + if(FULLTEST) { + ## special cases + nndist(X[FALSE]) + nndist(X[1]) + nndist(X[1:3], k=4) + nndist(X[1:3], k=1:4) + nnwhich(X[FALSE]) + nnwhich(X[1]) + nnwhich(X[1:3], k=4) + nnwhich(X[1:3], k=1:4) + nncross(X[1:3], Y[FALSE]) + nncross(X[1:3], Y[1]) + nncross(X[1:3], Y[1:3], k=4) + nncross(X[1:3], Y[1:3], k=1:4) + } + + ## ....... m dimensions ................ + + if(ALWAYS) { + rx <- function(n) { + B <- boxx(c(0,1),c(0,1),c(0,1),c(0,1)) + df <- replicate(4, runif(n), simplify=FALSE) + names(df) <- letters[23:26] + ppx(as.data.frame(df), B) + } + ## X <- runifpointx(42, B) + ## Y <- runifpointx(50, B) + X <- rx(42) + Y <- rx(50) + Y <- Y[nndist(Y) > 0.02] + DXY <- crossdist(X,Y) + + nn <- nndist(X) + nnP <- f(pairdist(X), 1) + if(any(abs(nn - nnP) > eps)) + stop("nndist.ppx does not agree with pairdist") + + nn5 <- nndist(X, k=5) + nn5P <- f(pairdist(X), 5) + if(any(abs(nn5 - nn5P) > eps)) + stop("nndist.ppx(k=5) does not agree with pairdist") + + nw <- nnwhich(X) + nwP <- g(pairdist(X), 1) + if(any(nw != nwP)) + stop("nnwhich.ppx does not agree with pairdist") + + nw5 <- nnwhich(X, k=5) + nw5P <- g(pairdist(X), 5) + if(any(nw5 != nw5P)) + stop("nnwhich.ppx(k=5) does not agree with pairdist") + + a <- nncross(X,Y) + ncd <- nncross(X,Y, what="dist") + ncw <- nncross(X,Y, what="which") + if(any(ncw != gg(DXY, 1))) + stop("incorrect result from nncross.ppx(what='which')") + a2 <- nncross(X,Y, k=2) + ncd <- nncross(X,Y, what="dist", k=2) + ncw <- nncross(X,Y, what="which", k=2) + if(any(ncw != gg(DXY, 2))) + stop("incorrect result from nncross.ppx(k=2, what='which')") + } + + if(FULLTEST) { + ## special cases + nndist(X[FALSE]) + nndist(X[1]) + nndist(X[1:3], k=4) + nndist(X[1:3], k=1:4) + nnwhich(X[FALSE]) + nnwhich(X[1]) + nnwhich(X[1:3], k=4) + nnwhich(X[1:3], k=1:4) + nncross(X[1:3], Y[FALSE]) + nncross(X[1:3], Y[1]) + nncross(X[1:3], Y[1:3], k=4) + nncross(X[1:3], Y[1:3], k=1:4) + } + + if(ALWAYS) { + ## test of agreement between nngrid.h and knngrid.h + ## dimyx=23 (found by trial-and-error) ensures that there are no ties + a <- as.matrix(nnmap(cells, what="which", dimyx=23)) + b <- as.matrix(nnmap(cells, what="which", dimyx=23, k=1:2)[[1]]) + if(any(a != b)) + stop("algorithms in nngrid.h and knngrid.h disagree") + + ## minnndist correctness + X <- redwood3 + eps <- sqrt(.Machine$double.eps) + mfast <- minnndist(X) + mslow <- min(nndist(X)) + if(abs(mfast-mslow) > eps) + stop("minnndist(X) disagrees with min(nndist(X))") + + ## maxnndist correctness + mfast <- maxnndist(X) + mslow <- max(nndist(X)) + if(abs(mfast-mslow) > eps) + stop("maxnndist(X) disagrees with max(nndist(X))") + } + + if(ALWAYS) { + ## minnndist, maxnndist code blocks + Y <- superimpose(amacrine, amacrine[10:20]) + a <- maxnndist(Y, positive=TRUE) + u <- maxnndist(Y, positive=TRUE, by=marks(Y)) + b <- minnndist(Y, positive=TRUE) + v <- minnndist(Y, positive=TRUE, by=marks(Y)) + + ## nnmap code blocks + A <- nnmap(cells[FALSE]) + A <- nnmap(cells, sortby="var") + A <- nnmap(cells, sortby="x") + A <- nnmap(cells, sortby="y") + B <- nnmap(cells[1:3], k=4) + B <- nnmap(cells[1:3], k=2:4) + D <- nnmap(cells, outputarray=TRUE) + } + + if(ALWAYS) { + #' tests for has.close() + #' (the default method uses nndist or pairdist, and can be trusted!) + a <- has.close(redwood, 0.05) + b <- has.close.default(redwood, 0.05) + if(any(a != b)) stop("Incorrect result for has.close(X, r)") + + a <- has.close(redwood, 0.05, periodic=TRUE) + a <- has.close.default(redwood, 0.05, periodic=TRUE) + if(any(a != b)) stop("Incorrect result for has.close(X, r, periodic=TRUE)") + + Y <- split(amacrine) + a <- with(Y, has.close(on, 0.05, off)) + b <- with(Y, has.close.default(on, 0.05, off)) + if(any(a != b)) stop("Incorrect result for has.close(X, r, Y)") + + a <- with(Y, has.close(on, 0.05, off, periodic=TRUE)) + b <- with(Y, has.close.default(on, 0.05, off, periodic=TRUE)) + if(any(a != b)) stop("Incorrect result for has.close(X, r, Y, periodic=TRUE)") + } + + if(ALWAYS) { + b <- bdist.pixels(letterR, style="coords") + d <- bdist.pixels(letterR, dimyx=64, method="interpreted") + } + + if(FULLTEST) { + ## nnfun.ppp + h <- nnfun(cells) + Z <- as.im(h) + d <- domain(h) + h <- nnfun(amacrine, value="mark") + d <- domain(h) + Z <- as.im(h) + h <- nnfun(longleaf, value="mark") + d <- domain(h) + Z <- as.im(h) + } + + }) Warning message: data contain duplicated points > > > proc.time() user system elapsed 1.25 0.15 1.40