R version 4.6.0 beta (2026-04-12 r89874 ucrt) -- "Because it was There" Copyright (C) 2026 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. > warn.tools_verbose <- FALSE > cat('Running test scripts:\n') Running test scripts: > library(levelSets) > pt_in_slice <- levelSets:::pt_in_slice > inRect <- levelSets:::inRect > flipRays <- levelSets:::flipRays > pt_on_line2 <- levelSets:::pt_on_line2 > warnBadResp <- levelSets:::warnBadResp > pairIdxs <- levelSets:::pairIdxs > linePosn0 <- levelSets:::linePosn0 > linePosn1 <- levelSets:::linePosn1 > lineVec <- levelSets:::lineVec > sliceMatch <- levelSets:::sliceMatch > internalOpt <- levelSets:::internalOpt > inpPts <- levelSets:::inpPts > inpPairs <- levelSets:::inpPairs > pairInfo <- levelSets:::pairInfo > respPts <- levelSets:::respPts > num_matrix <- levelSets:::num_matrix > naF <- levelSets:::naF > cat('axisRays.r', sep='\n') axisRays.r > # User-level functions to create standard sets of rays. > # > # axisRays > # randomRays > # flipRays (not public) > > oldopt <- lsetPkgOpt(reset=TRUE) > #----- > cat(" axisRays ... ") axisRays ... > #== Testing > > #----- 1D example. > spec <- fnSpec(inpnames="X1") > scl <- inpScale(matrix(5, nrow=1, ncol=1), spec=spec) > rays <- axisRays(spec, scale=scl) > stopifnot(dim(rayVec(rays)) == c(1, 1), rayVec(rays) == sqrt(5)) > # unit vector transformed to length sqrt(5) > > #----- 3D example to illustrate 'slice'ing. > #layout(matrix(c(1:4), nrow=2, byrow=TRUE)) > spec <- fnSpec(inpnames=c("X1", "X2", "X3")) > slice <- c(NA, NA, -2) > rays <- axisRays(spec, degree=2, slice=slice, origin=c(0, 0, -2)) # 4 rays > #print(names(rays)) > #plot(rays, eqAxes=TRUE) > stopifnot(all(pt_in_slice(rayPosn1(rays), slicemat=rbind(slice), sliceidx=1, + inptol=inpTol(spec))), + length(rays) == 4) > > # Specifying 'rect' and putting the ray origin on its boundary can lead > # to ray reversal or dropping rays: > rect <- inpRect(rbind(c(-1, 0, -2), c(1, 2, 3)), + refpt=c(1, 0, -2), # upper bdry of X1, lower of X2, X3 + spec=spec) > tools::assertWarning( + rays <- axisRays(spec, degree=2, slice=slice, rect=rect) # origin on bdry + # of 'rect' causes 'X1', 'X1(-)X2' rays to be reversed + # and 'X1(+)X2' ray to be dropped + , verbose=warn.tools_verbose) > #print(names(rays)) > stopifnot(all(pt_in_slice(rayPosn1(rays), slicemat=rbind(slice), sliceidx=1, + inptol=inpTol(spec)))) > stopifnot(all(inRect(rayPosn1(rays), rect=rect, bdryOnly=TRUE)), + length(rays) == 3) > #plot(rays, eqAxes=TRUE, main="Ray origin on boundary of 'rect'") > > # Can also specify 'scale', which is applied before ray vectors are > # extended to the boundary of 'rect'. > scl <- inpScale(rbind(c(1, 0.5, 0.5), c(0.5, 2, 0.5), c(0.5, 0.5, 4)), + spec=spec) # unequal scales, non-zero corr > rays <- axisRays(spec, degree=2, origin=c(1, 0, -2), slice=slice, scale=scl) > #plot(rays, eqAxes=TRUE, main="Scaled and rotated rays ...") > #plot(rect, add=TRUE) > # Rays 1 and 4 and their reverses are outside 'rect': > (chk <- inRect(rayPosn1(c(rays, flipRays(rays))), rect=rect)) rot1 rot2 rot1(+)rot2 rot1(-)rot2 rev_rot1 FALSE TRUE TRUE FALSE FALSE rev_rot2 rev_rot1(+)rot2 rev_rot1(-)rot2 FALSE FALSE FALSE > stopifnot(identical(unname(chk), c(FALSE, TRUE, TRUE, FALSE, + FALSE, FALSE, FALSE, FALSE))) > # So forcing boundary of 'rect' creates only 2 rays: > tools::assertWarning( + rays <- axisRays(spec, degree=2, slice=slice, rect=rect, scale=scl) + , verbose=warn.tools_verbose) > stopifnot(inherits(rays, "inpRays"), length(rays) == 2, + inRect(rayPosn1(rays), rect=rect, bdryOnly=TRUE)) > #plot(rays, eqAxes=TRUE, main="... confined to 'rect'") > # Shows how two rays are extended to bdry of 'rect' > > #layout(matrix(1)) > cat("OK\n") OK > > #----- > cat(" randomRays ... ") randomRays ... > #== Testing > spec <- fnSpec(inpnames=c("X")) # one-dimensional parameter space > set.seed(1) > rays <- randomRays(5, spec=spec) # random +/-1's > stopifnot(inherits(rays, "inpRays"), length(rays) == 5, + all(rayPosn1(rays) %in% c(-1, 1))) > posn1 <- rayPosn1(randomRays(5, spec=spec, origin=c(5))) # random 4's and 6's > stopifnot(is.matrix(posn1), dim(posn1) == c(5, 1), all(posn1 %in% c(4, 6))) > pRect <- inpRect(matrix(c(-1, 4), nrow=2), spec=spec) > posn1b <- rayPosn1(randomRays(5, rect=pRect)) # 'spec' taken from 'rect' > # random -1's and 4's > stopifnot(length(posn1b) == 5, all(posn1b %in% c(-1, 4))) > > #== 2-D examples. > #layout(matrix(1:4, nrow=2, byrow=TRUE)) > #----- Rays pass through random points on circle centered at (3,2) with > # radius 1: > spec <- fnSpec(inpnames=c("A", "B")) > set.seed(1) > rslt <- randomRays(50, spec=spec, origin=c(3, 2)) > #plot(rslt, eqAxes=TRUE, main="Random rays, unit length") > > # Effect of 'scale': > # Ellipsoid with wider range for A axis relative to B (semi-major and > # semi-minor axis lengths given by sqrt(diag(mat))). > scl <- inpScale(matrix(c(5, 0, 0, 1), nrow=2), spec=spec) > set.seed(1) > rslt <- randomRays(50, origin=c(3, 2), scale=scl) > #plot(rslt, eqAxes=TRUE, main="Ray lengths scaled to an ellipse") > > # Ellipsoid with wider range for A axis relative to B, and rotated: > scl <- inpScale(matrix(c(5, 2, 2, 1), nrow=2), spec=spec) > set.seed(1) > rslt <- randomRays(50, origin=c(3, 2), scale=scl) > #plot(rslt, eqAxes=TRUE, main="Rays scaled and rotated") > > # ... adjust lengths to put position 1.0 on the boundary of a rectangle: > rect <- inpRect(rbind(c(1, 1), c(5, 2.5)), spec=spec) > set.seed(1) > rslt <- randomRays(50, origin=c(3, 2), scale=scl, rect=rect) > #plot(rslt, eqAxes=TRUE, main="Ray position 1 on a 'rect' boundary") > > #-- When the sides of 'rect' have very different lengths, 'scale' makes > # a big difference. > #layout(matrix(c(1:3, 0), nrow=2, byrow=TRUE)) > # To highlight the effect, create a long, narrow rect, and use an off-center > # 'origin'. > corners <- rbind(c(0.5, 3), c(5.5, 3.5)) > rect <- inpRect(corners, spec=spec) > # Default for 'scale' adjusts for the side lengths, so rays fill 'rect': > set.seed(1) > rslt <- randomRays(50, origin=c(1, 3.25), rect=rect) > #plot(rslt, eqAxes=FALSE, main="Default 'scale' adjusts for 'rect' sides") > > # If we intentionally ignore differences in side lengths, rays may be sparse > # in the long directions: > set.seed(1) > rslt <- randomRays(50, origin=c(1, 3.25), rect=rect, + scale=inpScale(spec=spec)) > #plot(rslt, eqAxes=FALSE, main="No adjustment for 'rect' sides") > > #-- When 'origin' is itself on the boundary of 'rect', direction vectors > # that point outside 'rect' and lead to degenerate rays are instead reversed. > set.seed(1) > rslt <- randomRays(50, origin=c(2, 3.5), rect=rect) > #plot(rslt, eqAxes=FALSE, main="Ray origin on boundary of 'rect'") > > #layout(matrix(1)) > cat("OK\n") OK > > #----- > cat(" flipRays ... ") flipRays ... > #-- 1D example > # Set up an input space. > spec <- fnSpec(inpnames=c("a")) > rays <- inpRays(rbind(0, 3), origin=c(1), spec=spec) > flipped <- flipRays(rays) > stopifnot(inherits(flipped, "inpRays"), dim(rayPosn1(flipped)) == c(2, 1), + as.vector(rayPosn1(flipped)) == c(2, -1), + as.vector(rayOrigin(flipped)) == c(1)) > > # Origin on rectangle boundary: > pRect <- inpRect(rbind(-1, 2), spec=spec) > rays <- inpRays(c(0), origin=c(2), rect=pRect, spec=spec) > tools::assertWarning( + flipped <- flipRays(rays), verbose=warn.tools_verbose) # omit 1 ray > stopifnot(inherits(flipped, "inpRays"), length(flipped) == 0) > > #-- 3D examples > spec <- fnSpec(inpnames=c("A", "B", "C")) > rays <- inpRays(rbind(c(1, 1, 1), c(1, 2, 3)), origin=c(0, -1, 1), spec=spec) > rays2 <- flipRays(rays) > # Without a reference rectangle, flipRays() is its own inverse: > stopifnot(isTRUE(all.equal(flipRays(rays2), rays))) > > # Edge case: > rays2 <- flipRays(rays[0]) > stopifnot(inherits(rays2, "inpRays"), length(rays2) == 0) > > # Include a rect. > pRect <- inpRect(rbind(c(0,0,0), c(1,1,1)), spec=spec) > rays <- inpRays(vecs=rbind(c(-1,-1,0), c(0,0,1), c(0, -1, 1)), + origin=c(0.5, 0.4, 0.3), rect=pRect, spec=spec) > frays <- flipRays(rays) > # Flipped rays must have same origin as original, lie on the same lines, > # and have position 1 on the boundary of 'pRect'. > stopifnot(identical(rayOrigin(rays), rayOrigin(frays))) > stopifnot(all(pt_on_line2(rayPosn1(frays), lineidx=seq_len(length(rays)), + x0=rayOrigin(rays), x1=rayPosn1(rays), + inptol=inpTol(rays)))) > stopifnot(all(inRect(rayPosn1(frays), rect=pRect, bdryOnly=TRUE))) > > # Include a rect, and put the ray origin on its boundary: > pRect <- inpRect(rbind(c(0,0,0), c(1,1,1)), refpt=c(1,1,0.5), + spec=spec) > rays <- inpRays(vecs=rbind(c(-1,-1,0), c(0,0,1), c(0, -1, 1)), + rect=pRect, spec=spec) > # (First ray points inside the rectangle, second lies along the > # edge at the intersection of two faces, third lies in one of the > # faces containing the origin.) > # First and last rays are dropped because after flipping their only > # intersection with the rectangle is the origin: > tools::assertWarning( + frays <- flipRays(rays), verbose=warn.tools_verbose) > stopifnot(length(frays) == 1, rayPosn1(frays) == c(1, 1, 0)) > stopifnot(identical(rayOrigin(rays), rayOrigin(frays))) > stopifnot(all(pt_on_line2(rayPosn1(frays), lineidx=2, + x0=rayOrigin(rays), x1=rayPosn1(rays), + inptol=inpTol(rays)))) > cat("OK\n") OK > > #----- > lsetPkgOpt(oldopt) > cat('bdryFromRays.r', sep='\n') bdryFromRays.r > # 17 Feb 2025 / Richard Raubertas > # 24 Feb 2025 : 'fnObj' no longer has 'spec' arg. > # 12 Feb 2025 : Move 4D test to 'test_extras' directory. > > oldopt <- lsetPkgOpt(reset=TRUE) > > #----- > cat(" bdryFromRays ... ") bdryFromRays ... > #== Testing > #----- 1D example: > respf <- function(x, inclGrad) { + resp <- -1*(x^2 + 2*x + 1) + if (inclGrad) structure(resp, gradient= -1*(2*x + 2)) else resp } > respF <- fnObj(c("X1"), respfn=respf, hasgrad=TRUE) > pRays <- inpRays(matrix(c(5), ncol=1), origin=-2, spec=respF) > > rslt <- bdryFromRays(respF, rays=pRays, lsetthresh=-4) > stopifnot(inherits(rslt, "lsetSegs"), length(rslt) == 1, + isTRUE(all.equal(as.vector(ptCoord(rslt)), c(-3, 1), + tolerance=1e-6))) > summ <- summary(rslt) > stopifnot(summ$nlines == 1, colSums(summ$lineInfo)[1:4] == c(1, 0, 2, 0)) > > # Extract lines/rays from an 'lsetSegs' object. > stopifnot(identical(inpLines(rslt), pRays)) > > # No contour crossing between default positions -1, 0, and 1: > pRays <- inpRays(matrix(c(0), ncol=1), origin=-1, spec=respF) > #-- Only a censored segment is found if extra test positions are turned off: > rslt <- bdryFromRays(respF, rays=pRays, lsetthresh=-4, + control=list(initPosns2=numeric(0))) > grad <- attr(respInfo(rslt, inclGrad=TRUE), "gradient") > stopifnot(as.vector(ptCoord(rslt)) == c(-2, 0), + is.matrix(grad), dim(grad) == c(2, 1), + identical(dimnames(grad), list(NULL, "X1")), + as.vector(grad) == c(2, -2)) > summ <- summary(rslt) > stopifnot(summ$nlines == 1, colSums(summ$lineInfo)[1:4] == c(1, 0, 0, 2)) > #-- Both contour crossings found with default extra test positions: > rslt <- bdryFromRays(respF, rays=pRays, lsetthresh=-4) > stopifnot(inherits(rslt, "lsetSegs"), length(rslt) == 1, + isTRUE(all.equal(as.vector(ptCoord(rslt)), c(-3, 1), + tolerance=1e-6))) > summ <- summary(rslt) > stopifnot(summ$nlines == 1, colSums(summ$lineInfo)[1:4] == c(1, 0, 2, 0)) > grad <- attr(respInfo(rslt, inclGrad=TRUE), "gradient") > stopifnot(is.matrix(grad), dim(grad) == c(2, 1), + identical(dimnames(grad), list(NULL, "X1")), + as.vector(grad) == c(4, -4)) > > # Edge case: no rays. > rslt <- bdryFromRays(respF, rays=pRays[0], lsetthresh=-4) > stopifnot(inherits(rslt, "lsetSegs"), length(rslt) == 0) > summ <- summary(rslt) > stopifnot(summ$nlines == 0, colSums(summ$lineInfo)[1:4] == c(0, 0, 0, 0)) > > # Contour crossing exactly at an 'initPosn': > pRays <- inpRays(matrix(c(-3), ncol=1), origin=-2, spec=respF) > rslt <- bdryFromRays(respF, rays=pRays, lsetthresh=-4) > stopifnot(inherits(rslt, "lsetSegs"), length(rslt) == 1, + isTRUE(all.equal(as.vector(ptCoord(rslt)), c(1, -3), + tolerance=1e-6))) > > #----- 2D Multimodal example (two elliptical contours). > Ex <- dbl_ellipse_2dEx > respF <- Ex$fobj > fspec <- fnSpec(respF) > > # Small search region. > refregion <- inpRect(rbind(c(-8, -8), c(8, 8)), spec=fspec) > set.seed(1) > pRays <- randomRays(70, spec=fspec, origin=c(2, 1), rect=refregion) > > # Default 'initPosns2' is enough to find boundary points on both ellipses ... > rslt <- bdryFromRays(respF, rays=pRays, lsetthresh=-40) > summ <- summary(rslt) # 0 censored, 21 of 70 rays have more than one segment > stopifnot(summ$nlines == 70, colSums(summ$lineInfo)[1:4] == c(91, 0, 182, 0)) > # ... but 'initPosns' alone is not. > rslt <- bdryFromRays(respF, rays=pRays, lsetthresh=-40, + control=list(initPosns2=numeric(0))) > summ <- summary(rslt) # 0 censored, 21 of 70 rays have more than one segment > stopifnot(summ$nlines == 70, colSums(summ$lineInfo)[1:4] == c(70, 0, 100, 40)) > > # Test 'segBoundingRect': > rng <- apply(ptCoord(rslt), 2, range) > bdryrect <- segBoundingRect(rslt) # just includes points, so matches range: > stopifnot(isTRUE(all.equal(unname(rng), unname(rectCorners(bdryrect))))) > # Expand censored points: > bdryrect2 <- segBoundingRect(rslt, expand=c(1.0, 2.0)) > stopifnot(round(prod(rectSize(bdryrect))) == 229, + round(prod(rectSize(bdryrect2))) == 817) > > #-- Suppose we make poor choices of 'rect'. > # (a) Far from level set. > refregion2a <- inpRect(rbind(c(14, -2), c(15, -1)), spec=fspec) > set.seed(1) > pRays2a <- randomRays(70, spec=fspec, rect=refregion2a) > rslt2a <- bdryFromRays(respF, rays=pRays2a, lsetthresh=-40) # empty > summ <- summary(rslt2a) > stopifnot(summ$nlines == 70, colSums(summ$lineInfo)[1:4] == c(0, 0, 0, 0)) > > # (b) Small and contained entirely in the level set. > refregion2b <- inpRect(rbind(c(-1, -1), c(0, 0)), spec=fspec) > set.seed(1) > pRays2b <- randomRays(70, spec=fspec, rect=refregion2b) > rslt2b <- bdryFromRays(respF, rays=pRays2b, lsetthresh=-40) > summ <- summary(rslt2b) # all censored > stopifnot(summ$nlines == 70, colSums(summ$lineInfo)[1:4] == c(70, 0, 0, 140)) > #plot(rslt2b, rect=refregion2b) > > #-- Test feasible point checking. > Ex <- dbl_ellipse_X1bnd_2dEx > # Feasible region: X1 >= 0. > respF <- Ex$fobj > rslt <- bdryFromRays(respF, rays=pRays, lsetthresh=-40) > summ <- summary(rslt) # 0 of 182 censored, 50 on feas bdry; 21 rays w/ multiple segs > stopifnot(summ$nlines == 70, colSums(summ$lineInfo)[1:4] == c(91, 50, 132, 0)) > > #-- Use a 'rayProfiles' object as 'lines' argument. > profs <- respProfiles(respF, lines=pRays, lsetthresh=-40, + posns=c(-2, -1.5, -1, -0.5, 0, 0.5, 1, 1.5, 2)) > #-- Just extract existing bdry intersections from 'profs': > rslt2 <- lsetSegs(profs) # '.respProfiles' method > summ <- summary(rslt2) > stopifnot(summ$nlines == 70, colSums(summ$lineInfo)[1:4] == c(102, 50, 148, 6)) > # One segment crosses the gap between ellipses and so is invalid. Can be > # fixed with 'lsetSegsCheck': > chk <- lsetSegsCheck(rslt2, fnobj=respF, posns=(1:4)/5) > stopifnot(length(chk$validIn) == length(rslt2), + sum(chk$validIn) == length(rslt2) - 1, + length(chk$validOut) == length(rslt2) + 1, all(chk$validOut)) > rslt2b <- chk$segs > summ <- summary(rslt2b) > stopifnot(summ$nlines == 70, colSums(summ$lineInfo)[1:4] == c(103, 50, 150, 6)) > > #-- Try changing 'lsetthresh' value from what was used in profiling. > rslt3b <- bdryFromRays(respF, rays=profs, lsetthresh=-30) > stopifnot(lsetThresh(rslt3b) == -30) > # Should be the same as if we started with a profile using the new threshold: > profs_alt <- respProfiles(respF, lines=pRays, lsetthresh=-30, + posns=c(-2, -1.5, -1, -0.5, 0, 0.5, 1, 1.5, 2)) > rslt3b_alt <- bdryFromRays(respF, rays=profs_alt) > stopifnot(isTRUE(all.equal(rslt3b, rslt3b_alt, tolerance=1e-6)), + isTRUE(all.equal(summary(rslt3b), summary(rslt3b_alt), + tolerance=1e-6))) > > #-- Ray profiles that never intersect the level set (e.g., any level set > # defined by a positive contour value): > profs <- respProfiles(respF, lines=pRays, posns=c(-1, 0, 1), lsetthresh=10) > rslt4 <- bdryFromRays(respF, rays=profs) # should be empty > stopifnot(inherits(rslt4, "lsetSegs"), length(rslt4) == 0) > summ <- summary(rslt4) > stopifnot(summ$nlines == 70, colSums(summ$lineInfo)[1:4] == c(0, 0, 0, 0)) > > # Test 'segBoundingRect' with no level set points: > tools::assertWarning({ bdryrect <- segBoundingRect(rslt4) } + ) > stopifnot(is.null(bdryrect)) > > cat("OK\n") OK > > #----- > lsetPkgOpt(oldopt) > > cat('bdrySearch.r', sep='\n') bdrySearch.r > # 24 Feb 2025 / Richard Raubertas > # 24 Jul 2025 : Update names of 'control' parameters. > # 27 Aug 2025 : Return component 'origins' renamed to 'srchOrigins'. > # 12 Feb 2026 : Remove 4D tests with t-distribution example (now in > # 'tests_extra'). > # 01 Mar 2026 : No longer reference 'lvlSets_env' for 2D example. > > oldopt <- lsetPkgOpt(reset=TRUE) > lsetPkgOpt(warn.showBadRespPts=0) > > #----- > cat(" bdrySearch ... \n") bdrySearch ... > #== Testing > > #----- 1D example. > cat(" 1D examples ... ") 1D examples ... > respf <- function(x, ...) { x*cos(x) } > #xx <- seq(-3, 10, by=0.25) > #plot(xx, respf(xx), type="l") > #abline(h=0.2, lty=2) > > fobj <- fnObj(c("X1"), respfn=respf) > rect1 <- inpRect(cbind(c(-3, 10)), spec=fobj) > srch1 <- bdrySearch(fobj, lsetthresh=0.2, rect=rect1) > bdry1 <- srch1$segs > summ <- summary(bdry1) > stopifnot(inherits(bdry1, "lsetSegs"), length(bdry1) == 2, # 2 segments + nrow(srch1$srchOrigins) == 1, # matrix not 'respPts' + summ$nlines == 1, colSums(summ$lineInfo)[1:4] == c(2, 0, 4, 0)) > # Upper segment is correct. Lower seg should really be > # many segs; multiple oscillations below X1=1 are missed: > if (FALSE) { + xx <- seq(-30, 10, by=0.25) + plot(xx, respf(xx), type="l") + abline(h=0.2, lty=2) + abline(v=as.vector(ptCoord(bdry1)), col="red") + } > > #-- Let search go wild chasing oscillations: allow multiple search steps > # even in 1D. > tools::assertWarning( + srch1b <- bdrySearch(fobj, lsetthresh=0.2, rect=rect1, + control=list(maxSteps=5)) + ) # 12 seg pairs reported inconsistent > # With 'initPosns2' values > 1 and 'rect' updating the search region expands > # dramatically over steps, and futilely, since the number of segments is > # infinite. > bdry1b <- srch1b$segs > summ <- summary(bdry1b) # 7 segs on one line, no censoring > stopifnot(inherits(bdry1b, "lsetSegs"), length(bdry1b) == 7, # 7 segments + nrow(srch1b$srchOrigins) == 5, # matrix not 'respPts' + summ$nlines == 1, colSums(summ$lineInfo)[1:4] == c(7, 0, 14, 0)) > #??????? Add option in lsetSegsCheck to only check validity or only consistency > temp <- lsetSegsCheck(bdry1b, fnobj=fobj, posns=(1:4)/5) # default 1 fix cycle > stopifnot(length(temp$validIn) == length(bdry1b), sum(temp$validIn) == 1, + length(temp$validOut) == 16, sum(temp$validOut) == 5, + is.matrix(temp$inconOut), nrow(temp$inconOut) == 20) > # Goes from 6/7 invalid to 11/16, with 20 inconsistent pairs > temp2 <- lsetSegsCheck(temp$segs, fnobj=fobj, posns=(1:4)/5, fixCycles=5) > stopifnot(length(temp2$validIn) == length(temp$segs), sum(temp2$validIn) == 5, + length(temp2$validOut) == 20, sum(temp2$validOut) == 20, + is.matrix(temp2$inconOut), nrow(temp2$inconOut) == 3) > # Goes from 11/16 invalid to 0/20, with 3 inconsistent pairs > # However, ~60 segs are expected given the range of seg box > > #-- Limit rect expansion by restricting 'initPosns2' to values <= 1. > tools::assertWarning( + srch1d <- bdrySearch(fobj, lsetthresh=0.2, rect=rect1, + control=list(maxSteps=5, initPosns2=c(0.25, 0.5))) + ) # 3 seg pairs reported inconsistent > bdry1d <- srch1d$segs > summ <- summary(bdry1d) > stopifnot(inherits(bdry1d, "lsetSegs"), length(bdry1d) == 3, # 3 segments + nrow(srch1d$srchOrigins) == 5, # matrix not 'respPts' + summ$nlines == 1, colSums(summ$lineInfo)[1:4] == c(3, 0, 6, 0)) > # Three (inconsistent) segs in a box only slightly larger than 'rect1' (due > # to multiple search origins). Correct answer is three segs, but two > # intermediate crossing pts are missed. > temp <- lsetSegsCheck(bdry1d, posn=(1:4)/5, fnobj=fobj, fixCycles=3) > stopifnot(length(temp$validIn) == length(bdry1d), sum(temp$validIn) == 1, + length(temp$validOut) == 3, sum(temp$validOut) == 3, + is.matrix(temp$inconOut), nrow(temp$inconOut) == 0) > #table(temp$validIn) # 2 of 3 orignal segs are invalid > #table(temp$validOut) # now 3 valid segs > #print(temp$inconOut) # none are inconsistent > #print(segInfo(temp$segs)) > if (FALSE) { + xx <- seq(-5, 10, by=0.25) + plot(xx, respf(xx), type="l") + abline(h=0.2, lty=2) + abline(v=as.vector(ptCoord(temp$segs)), col="red") + } > # All segs correct. > > cat("OK\n") OK > > #----- 2D multimodal example (level set has two separate elliptical parts), > cat(" 2D examples ... ") 2D examples ... > Ex <- dbl_ellipse_X1bnd_2dEx > # Feasible region: X1 >= 0. > respF <- Ex$fobj > fspec <- fnSpec(respF) > > # Initial search region (default reference point is center, (0, 0)): > refregion1 <- inpRect(rbind(c(-8, -8), c(8, 8)), spec=fspec) > > # Boundary search with defaults: > srch1 <- bdrySearch(respF, lsetthresh=-40, rect=refregion1) > bdry1 <- srch1$segs > summ <- summary(bdry1) > stopifnot(inherits(bdry1, "lsetSegs"), length(bdry1) == 9, + nrow(srch1$srchOrigins) == 4, # matrix not 'respPts' + summ$nlines == 8, colSums(summ$lineInfo)[1:4] == c(9, 3, 15, 0)) > > # Use the 'control' argument to increase the number of segments. > srch1b <- bdrySearch(respF, lsetthresh=-40, rect=refregion1, + control=list(axisDegree=2, # increases rays per origin + maxSteps=8)) > bdry1b <- srch1b$segs > srch1b$srchOrigins # 8 origins used, 4 rays per origin X1 X2 [1,] 0.000000 0.000000 [2,] 6.611926 -1.938269 [3,] 7.542630 8.764301 [4,] 13.534678 5.643419 [5,] 2.518875 -4.608141 [6,] 3.383670 10.930442 [7,] 8.796400 4.946464 [8,] 3.084459 1.686009 > summ <- summary(bdry1b) # 34 segments on 32 lines; 0 censored endpts; 2 lines w/ > 1 seg > stopifnot(inherits(bdry1b, "lsetSegs"), length(bdry1b) == 34, + nrow(srch1b$srchOrigins) == 8, # matrix not 'respPts' + summ$nlines == 32, colSums(summ$lineInfo)[1:4] == c(34, 15, 53, 0)) > #plot(bdry1b) # many segs cross the gap between ellipses (i.e., invalid) > > # Check segment validity, and attempt to fix any invalid ones. > chk <- lsetSegsCheck(bdry1b, fnobj=respF, posns=(1:4)/5) > stopifnot(table(chk$validIn) == c(11, 23), # 11 invalid segments + table(chk$validOut) == c(45), # Fixed by splitting each into two + is.matrix(chk$inconOut), nrow(chk$inconOut) == 0) > bdry1c <- chk$segs > summ <- summary(bdry1c) > stopifnot(inherits(bdry1c, "lsetSegs"), length(bdry1c) == 45, + summ$nlines == 32, colSums(summ$lineInfo)[1:4] == c(45, 15, 75, 0)) > > # Fewer invalid segments if we use a finer and wider grid of starting > # positions for search. > srch3 <- bdrySearch(respF, lsetthresh=-40, rect=refregion1, + currentBdry=srch1, + control=list(axisDegree=2, initPosns=(0:8)/4)) > stopifnot(nrow(srch3$srchOrigins) == 8) > bdry3 <- srch3$segs > chk <- lsetSegsCheck(bdry3, fnobj=respF, posns=(1:4)/5) > stopifnot(table(chk$validIn) == c(4, 26), # Four invalid segments + table(chk$validOut) == c(34)) # Fixed by splitting each into two > bdry3b <- chk$segs > summ <- summary(bdry3b) > stopifnot(inherits(bdry3b, "lsetSegs"), length(bdry3b) == 34, + summ$nlines == 24, colSums(summ$lineInfo)[1:4] == c(34, 9, 59, 0)) > > # Try skipping the random rotation of ray vectors. > srch4 <- bdrySearch(respF, lsetthresh=-40, rect=refregion1, + control=list(axisDegree=2, + maxSteps=8, rotateRays=FALSE)) > bdry4 <- srch4$segs > summ <- summary(bdry4) > stopifnot(inherits(bdry4, "lsetSegs"), length(bdry4) == 35, + nrow(srch4$srchOrigins) == 8, # matrix not 'respPts' + summ$nlines == 32, colSums(summ$lineInfo)[1:4] == c(35, 12, 58, 0)) > #plot(bdry4, rect=refregion1) > # Not very good. > chk <- lsetSegsCheck(bdry4, fnobj=respF, posns=(1:4)/5) > stopifnot(table(chk$validIn) == c(7, 28), table(chk$validOut) == 42, + nrow(chk$inconOut) == 0) > #plot(chk$segs, rect=refregion1) # not bad > > # Conclusion: > # a) Important to have right search region ('rect'), updating after initial > # guess(es) if necessary. > # b) Need grid of 'initPosns' fine enough to match scale of interesting > # features of level set shape (concavity or multiple parts). > # c) Need to increase number of boundary pairs to see detailed shape of > # level set. Can be done by increasing the number of origins > # ('maxSteps') initially or over several calls, and/or by > # increasing 'axisDegree'. > # d) 'rotateRays' seems to be useful. > > cat("OK\n") OK > > #----- > lsetPkgOpt(oldopt) > cat('fnObj.r', sep='\n') fnObj.r > # fnObj > # hasGrad > # warnBadResp* > # fnObjCheck* (???? incomplete and disabled) > > # *= not public > > # 12 Feb 2025 / Richard Raubertas > # 24 Feb 2025 : 'fnObj' no longer has a 'spec' arg. Test new 'update.fnObj'. > > oldopt <- lsetPkgOpt(reset=TRUE) > lsetPkgOpt(warn.showBadRespPts=0) # avoid printing matrices, data frames > #----- > cat(" fnObj ... ") fnObj ... > #== Testing > # Function of 3 input variables, with gradient available, and an > # extra argument 'powr'. > f <- function(x, inclGrad=FALSE, powr) { + y <- x[, "X1"]^powr[1] + x[, "X2"]^powr[2] + x[, "X3"]^powr[3] + if (inclGrad) { + grad <- array(NA_real_, dim=dim(x), dimnames=dimnames(x)) + for (j in 1:3) { + grad[, j] <- { if (powr[j] == 0) 0 + else powr[j] * (x[, j]^(powr[j] - 1)) } + } + attr(y, "gradient") <- grad + } + y + } > # Associated feasibility function, passed the same extra argument. > feasf <- function(x, powr) { + feas <- rep(TRUE, nrow(x)) + # Want to avoid applying non-integer powers to negative numbers. + fracpow <- (powr != round(powr)) + for (j in which(fracpow)) { + feas <- feas & (x[, j] >= 0) # NOTE: '>=', not '>'. Feasible set + # should be _closed_. + } + feas + } > > fnobj <- fnObj(c("X1", "X2", "X3"), respfn=f, feasfn=feasf, + hasgrad=TRUE, powr=c(-1, 0, 2.5)) > stopifnot(hasGrad(fnobj)) > fspec <- fnSpec(fnobj) > > # Tests: pts with NA coord, or generating Inf resp: > pts <- rbind(c(0, 0, 0), c(2, 0, -1), c(1, 6, 1), c(NA, 1, 1)) > colnames(pts) <- inpNames(fnobj) > tools::assertError(respInfo(pts, fnobj=fnobj), # NA pt coord not allowed + verbose=warn.tools_verbose) > pts[4, ] <- c(-1, 1, 1) > tools::assertWarning(respInfo(pts, fnobj=fnobj), # Inf response + verbose=warn.tools_verbose) > > pts[1, ] <- c(3, -2, 0) > resp1 <- respInfo(pts, fnobj=fnobj, inclGrad=TRUE) > stopifnot(is.data.frame(resp1), nrow(resp1) == nrow(pts), + resp1$feas == c(TRUE, FALSE, TRUE, TRUE), + is.na(resp1$resp) == !resp1$feas) > grad <- attr(resp1, "gradient") > stopifnot(dim(grad) == dim(pts), + apply(grad, 1, function(y) {any(is.na(y))}) == !resp1$feas) > resp2 <- respInfo(pts, fnobj=fnobj, inclGrad=FALSE) > grad <- attr(resp2, "gradient") > stopifnot(is.null(grad)) > > # Test '.fnSpec' method: > fspec2 <- fnSpec(inpnames=c("X1", "X2", "X3")) > fnobj2 <- fnObj(fspec2, respfn=f, feasfn=feasf, + hasgrad=TRUE, powr=c(-1, 0, 2.5)) > resp2 <- respInfo(pts, fnobj=fnobj2, inclGrad=TRUE) > stopifnot(identical(fnSpec(fnobj), fnSpec(fnobj2)), + identical(resp1, resp2)) > > #-- Feasibility function based on 'feasbnds'. > fnobj <- fnObj(c("X1", "X2", "X3"), respfn=f, + feasbnds=rbind(rep(0, 3), rep(NA, 3)), hasgrad=FALSE, + powr=c(-1, 0, 2.5), warnInfresp=FALSE) > stopifnot(!hasGrad(fnobj)) > pts <- rbind(c(0, 0, 0), c(2, 0, -1), c(1, 6, 1), c(-2, 1, 1)) > colnames(pts) <- inpNames(fnobj) > resp3 <- respInfo(pts, fnobj=fnobj) # no warning about Inf response > stopifnot(is.data.frame(resp3), nrow(resp3) == nrow(pts), + resp3$feas == c(TRUE, FALSE, TRUE, FALSE), + is.na(resp3$resp) == !resp3$feas) > > #-- Combine explicit feasibility function with 'feasbnds'. > fnobj <- fnObj(c("X1", "X2", "X3"), respfn=f, feasfn=feasf, + feasbnds=list(X2=c(2, NA)), hasgrad=TRUE, powr=c(-1, 0, 2.5), + warnInfresp=0) > pts <- rbind(c(0, 0, 0), c(2, 3, -1), c(1, 6, 1), c(0, 4, 1), c(1, 1, 1)) > colnames(pts) <- inpNames(fnobj) > # Negative value makes 2nd pt infeasible, lower bound makes 1st and 5th > # infeasible. > resp4 <- respInfo(pts, fnobj=fnobj) # no warning about Inf response > stopifnot(is.data.frame(resp4), nrow(resp4) == nrow(pts), + resp4$feas == c(FALSE, FALSE, TRUE, TRUE, FALSE), + is.na(resp4$resp) == !resp4$feas) > grad <- attr(resp4, "gradient") > stopifnot(is.null(grad)) # b/c default 'inclGrad' is FALSE even w/ 'hasGrad=TRUE' > > #----- Testing 'warnNAresp', and 'update' method. > respf <- function(x, inclGrad=FALSE) { + resp <- rowSums(x) + ifelse(x[, 1] + x[, 2] >= 2, resp, NA) # implicit feasible region + } > feasf <- function(x) { # explicit feasible region + (x[, 1] - x[, 2] >= 1) + } > pts <- rbind(c(3, 1, 4), # feasible by all tests + c(5, -1, 4), # not feas wrt box constraint + c(1, 1, 4), # not feas wrt explicit constraint + c(1, 0, 4)) # not feas wrt implicit constraint > # No feasfn, silently accept implicitly infeas pts: > fnobj <- fnObj(c("X1", "X2", "X3"), respfn=respf, warnNAresp=FALSE, + feasbnds=list(X2=c(0, NA))) # box constraint feasible region > resp <- respInfo(pts, fnobj=fnobj) > stopifnot(identical(resp[, "resp"], c(8, NA, 6, NA)), + resp[, "feas"] == !is.na(resp$resp), + identical(resp[, "lset"], c(NA, FALSE, NA, FALSE))) > # No feasfn, warn about implicitly infeas pts: > fnobj2 <- fnObj(c("X1", "X2", "X3"), respfn=respf, warnNAresp=TRUE, + feasbnds=list(X2=c(0, NA))) > tools::assertWarning(resp2 <- respInfo(pts, fnobj=fnobj2), + verbose=warn.tools_verbose) > stopifnot(identical(resp2[, "resp"], c(8, NA, 6, NA)), + resp2[, "feas"] == !is.na(resp2$resp), + identical(resp2[, "lset"], c(NA, FALSE, NA, FALSE))) > # No feasfn, treat implicitly infeas pts as error: > fnobj <- update(fnobj, warnNAresp=2) > tools::assertError(resp <- respInfo(pts, fnobj=fnobj), + verbose=warn.tools_verbose) > # Same, relying on default 'warnNAresp': > fnobj <- fnObj(c("X1", "X2", "X3"), respfn=respf) > tools::assertError(resp <- respInfo(pts, fnobj=fnobj), + verbose=warn.tools_verbose) > > # Specify feasfn, silently accept implicitly infeas pts: > fnobj <- fnObj(c("X1", "X2", "X3"), respfn=respf, feasfn=feasf, + feasbnds=list(X2=c(0, NA)), warnNAresp=FALSE) > resp <- respInfo(pts, fnobj=fnobj) > stopifnot(identical(resp[, "resp"], c(8, NA, NA, NA)), + resp[, "feas"] == !is.na(resp$resp), + identical(resp[, "lset"], c(NA, FALSE, FALSE, FALSE))) > # Specify feasfn, warn about implicitly infeas pts: > fnobj <- update(fnobj, warnNAresp=TRUE) > tools::assertWarning(resp <- respInfo(pts, fnobj=fnobj), + verbose=warn.tools_verbose) > stopifnot(identical(resp[, "resp"], c(8, NA, NA, NA)), + resp[, "feas"] == !is.na(resp$resp), + identical(resp[, "lset"], c(NA, FALSE, FALSE, FALSE)), + fnobj$warnNAresp == 1) > # Specify feasfn, treat implicitly infeas pts as error: > fnobj <- update(fnobj, warnNAresp=2) > tools::assertError(resp <- respInfo(pts, fnobj=fnobj), + verbose=warn.tools_verbose) > > # No explicit feasible region at all, just implicit one. > respf <- function(x) { + resp <- rowSums(x) + ifelse(x[, 1] + x[, 2] >= 2, resp, NA) # implicit feasible region + } > # Default issues error: > fnobj <- fnObj(c("X1", "X2", "X3"), respfn=respf) > tools::assertError(resp <- respInfo(pts, fnobj=fnobj), + verbose=warn.tools_verbose) > # User must intentionally disable error to rely on implicit feasible > # region: > fnobj <- update(fnobj, warnNAresp=0) > resp <- respInfo(pts, fnobj=fnobj) > stopifnot(identical(resp[, "resp"], c(8, 8, 6, NA)), + resp[, "feas"] == !is.na(resp$resp), + "feasFn" %in% names(fnobj) && is.null(fnobj$feasFn)) > stopifnot(!hasGrad(fnobj)) # default > > cat("OK\n") OK > > #----- > cat(" warnBadResp ... ") warnBadResp ... > #== Testing > pts <- rbind(c(1, 2), c(2, 3), c(3, 4), c(4, 5), c(5, 6)) > resp <- c(4, Inf, -Inf, 2, NA) > tools::assertWarning(warnBadResp(resp, pts, type="Inf"), + verbose=warn.tools_verbose) # 2 pts printed > tools::assertWarning(warnBadResp(resp, pts, type="NA"), + verbose=warn.tools_verbose) # 1 pt printed > #lsetPkgOpt(warn.showBadRespPts=1) > tools::assertWarning(warnBadResp(resp, pts, type="Inf"), + verbose=warn.tools_verbose) # 1 pt printed > #lsetPkgOpt(warn.showBadRespPts=0) > tools::assertWarning(warnBadResp(resp, pts, type="Inf"), + verbose=warn.tools_verbose) # no pts printed > > #lsetPkgOpt(warn.showBadRespPts=1) > tools::assertError(warnBadResp(resp, pts, type="Inf", level=2), + verbose=warn.tools_verbose) # 1 pt printed > > cat("OK\n") OK > > #----- > cat(" update.fnobj ... ") update.fnobj ... > fnobj <- fnObj(c("X1", "X2"), respfn=function(x) { colSums(x) }, + inptol=c(1e-6, 1e-8)) > inptol1 <- inpTol(fnobj) > stopifnot(is.matrix(inptol1), dim(inptol1) == c(2, 2), + unname(inptol1) == cbind(c(1e-6, 1e-8), c(1e-6, 1e-8))) > > fnobj2 <- update(fnobj, inptol=c(1e-4, 2e-8)) > inptol2 <- inpTol(fnobj2) > stopifnot(is.matrix(inptol2), dim(inptol2) == c(2, 2), + unname(inptol2) == cbind(c(1e-4, 2e-8), c(1e-4, 2e-8))) > > cat("OK\n") OK > > #----- > lsetPkgOpt(oldopt) # restore original options > cat('fnSpec.r', sep='\n') fnSpec.r > # Define 'fnSpec' class and related accessor functions. > > # fnSpec > # print.fnSpec* > # fnspec_accessor* > # inpDim > # inpNames > # feasBnds > # inpTol > # respTol > > # *= not public > > # 11 Feb 2025 / Richard Raubertas > > oldopt <- lsetPkgOpt(reset=TRUE) > #----- > cat(" fnSpec ... ") fnSpec ... > #== Testing > fspec <- fnSpec(inpnames=c("X1", "X2", "X3"), + feasbnds=list(X2=c(0, 10), X3=c(NA, 5))) > stopifnot(identical(fspec, fnSpec(fspec))) > stopifnot(dim(inpTol(fspec)) == c(2, 3)) > stopifnot(!is.matrix(respTol(fspec)), length(respTol(fspec)) == 2) > stopifnot(inpDim(fspec) == 3) > stopifnot(inpNames(fspec) == c("X1", "X2", "X3")) > stopifnot(feasBnds(fspec) == rbind(c(-Inf, 0, -Inf), c(Inf, 10, 5))) > > # Updating: > fspec <- fnSpec(inpnames=c("X1", "X2", "X3")) > fspec2 <- fnSpec(fspec, feasbnds=list(X2=c(0, 10), X1=c(NA, 5))) > stopifnot(identical(inpTol(fspec2), inpTol(fspec))) > fspec2b <- fnSpec(fspec, feasbnds=feasBnds(fspec2)) > stopifnot(identical(fspec2, fspec2b)) > fspec3 <- fnSpec(fnSpec(fspec2, resptol=c(1e-6, 1e-7)), + resptol=respTol(fspec2)) > stopifnot(identical(fspec3, fspec2)) > > cat("OK\n") OK > > #----- > lsetPkgOpt(oldopt) > cat('inpLines.r', sep='\n') inpLines.r > # inpLines (S3 generic) > # inpLines.default (handles matrices of pt coords, and 'inpPts' objects) > # inpLines.inpPairs* (also handles objects inheriting from 'inpPairs') > # inpLines.inpLines* (also handles objects inheriting from 'inpLines') > # inpLines.fnSpec (allows creation of empty 'inpLines' object) > # inpLines.respProfiles > # inpLines.lsetSegs (needed, otherwise 'inpPairs' method takes precedence > # over 'respProfiles' method) > > # *= not public > > oldopt <- lsetPkgOpt(reset=TRUE) > #----- > cat(" inpLines ... ") inpLines ... > #== Testing > > fspec <- fnSpec(inpnames=c("X1", "X2", "X3")) > pts <- rbind(c(0, 1, 3), c(2, 2, 4), c(5, -1, 2), c(1, 1, 1)) > # The default method of defining positions 0 and 1 for each line is > # to use consecutive pairs of rows of 'pts': > ilines <- inpLines(pts, spec=fspec) > stopifnot(length(ilines) == 2) > > # Alternative way of specifying positions 0 and 1 is by two matrices: > ilines2 <- inpLines(pts[c(2, 4), ], x0=pts[c(1, 3), ], spec=fspec) > stopifnot(identical(ptCoord(ilines), ptCoord(ilines2, which="paired"))) > > names(ilines) <- c("foo", "bar") > stopifnot(colnames(pairIdxs(ilines)) == c("foo", "bar")) > # Downshifting to 'inpPairs' should remove names. > stopifnot(is.null(colnames(pairIdxs(inpPairs(ilines))))) > > # Empty 'inpLines' object. > ilines0 <- inpLines(fspec) > stopifnot(identical(class(ilines0), class(ilines)), + length(ilines0) == 0, + identical(ilines0, ilines[0])) > ilines00 <- inpLines(pts[0, , drop=FALSE], spec=fspec) > stopifnot(identical(ilines0, ilines00)) > > # Non-default pairing info. > pidxs <- rbind(c(1, 1, 1), c(2, 3, 4)) > colnames(pidxs) <- LETTERS[1:3] > ilines2 <- inpLines(pts, pairidxs=pidxs, spec=fspec) > # Accessor functions: > stopifnot(length(ilines2) == 3) > stopifnot(dim(pairIdxs(ilines2)) == c(2, 3), + colnames(pairIdxs(ilines2)) == LETTERS[1:3]) > pairInfo(ilines2) # data frame pairType idx1 idx2 1 1 2 2 1 3 3 1 4 > stopifnot(identical(unname(linePosn0(ilines2)), pts[c(1, 1, 1), ])) > stopifnot(identical(unname(linePosn1(ilines2)), pts[2:4, ])) > stopifnot(isTRUE(all.equal(lineVec(ilines2), + linePosn1(ilines2) - linePosn0(ilines2)))) > stopifnot(dim(ptCoord(ilines2)) == c(nrow(pts), inpDim(fspec))) > # 4 rows (raw point coords, ignoring pairing) > stopifnot(dim(ptCoord(ilines2, which="paired")) == + c(2*length(ilines2), inpDim(fspec))) # 2*number of pairs rows > names(ilines2) <- c("foo", "bar", "baz") > stopifnot(identical(linePosn0(ilines2), ptCoord(ilines2, which="1st"))) > stopifnot(identical(linePosn1(ilines2), ptCoord(ilines2, which="2nd"))) > > # Combine multiple 'inpLines' objects. > stopifnot(identical(ilines2, c(ilines2))) > # Internal representation will not be identical, so this may be FALSE ... > #identical(ilines2, c(ilines2[1:2], ilines2[0], ilines2[3], NULL)) > # ... but actual point pairs will be the same: > ptprs <- ptCoord(ilines2, which="paired") > cptprs <- ptCoord(c(ilines2[1:2], ilines2[0], ilines2[3], NULL), + which="paired") > stopifnot(identical(ptprs, cptprs)) > > infoprs <- pairInfo(ilines2) > cinfoprs <- pairInfo(c(ilines2[1:2], ilines2[0], ilines2[3], NULL)) > stopifnot(identical(infoprs$foo, cinfoprs$foo)) > > # When 'x' is already an 'inpLines' object. > stopifnot(identical(ilines, inpLines(ilines)), + identical(ilines0, inpLines(ilines0))) > > # 'strictClass' argument. > rays <- inpRays(pts, origin=c(2, 2, 2), spec=fspec) > lines_from_rays <- inpLines(rays) # default 'strictClass=TRUE' > stopifnot(!inherits(lines_from_rays, "inpRays")) > still_rays <- inpLines(rays, strictClass=FALSE) > stopifnot(identical(rays, still_rays), + identical(ilines, inpLines(ilines, strictClass=FALSE))) > > # Degenerate lines. > pts <- rbind(c(0, 1, 3), c(2, 2, 4), c(5, -1, 2), c(0, 1, 3) + 1e-8) > # One pair tied exactly, another to within default tolerance: > pairs <- rbind(c(1, 1, 1, 1), c(1, 2, 3, 4)) > # Default is to drop tied pairs and warn: > tools::assertWarning(lines <- inpLines(pts, pairidxs=pairs, spec=fspec), + verbose=warn.tools_verbose) > # (Only the 1st warning message is printed, but both tied pairs are dropped:) > stopifnot(inherits(lines, "inpLines"), length(lines) == 2) > # Drop tied pairs silently: > lines <- inpLines(pts, pairidxs=pairs, spec=fspec, warnDegenerate=0) > stopifnot(inherits(lines, "inpLines"), length(lines) == 2) > # Make tied pairs an error: > tools::assertError(lines <- inpLines(pts, pairidxs=pairs, spec=fspec, + warnDegenerate=2), + verbose=warn.tools_verbose) > cat("OK\n") OK > > #----- > lsetPkgOpt(oldopt) > cat('inpObjAccessors.r', sep='\n') inpObjAccessors.r > > # rectCorners > # rectRefPt > # rectSize > # sliceMatch (not public) > > oldopt <- lsetPkgOpt(reset=TRUE) > #----- > cat(" sliceMatch ... ") sliceMatch ... > #== Testing > fspec <- fnSpec(inpnames=c("a", "b")) > inptol <- inpTol(fspec) > slices <- sliceMat(rbind(c(1, NA), c(NA, 4)), spec=fspec) > y <- rbind(c(1, 3), c(2, 2), c(1, 4), c(3, 4)) > stopifnot(identical(sliceMatch(slices, y, inptol=inptol), + list(`a=1`=c(1L, 3L), `b=4`=c(3L, 4L)))) > stopifnot(identical(sliceMatch(slices, inpPts(y, spec=fspec), inptol=inptol), + list(`a=1`=c(1L, 3L), `b=4`=c(3L, 4L)))) # same > stopifnot(identical(sliceMatch(slices, inpPairs(y, spec=fspec), inptol=inptol), + list(`a=1`=integer(0), `b=4`=c(2L)))) > stopifnot(identical(sliceMatch(slices, inpLines(y, spec=fspec), inptol=inptol), + list(`a=1`=integer(0), `b=4`=c(2L)))) # same > > cat("OK\n") OK > > #----- > lsetPkgOpt(oldopt) > > cat('inpObjClasses.r', sep='\n') inpObjClasses.r > # Classes of objects in an input space. ('inpLines' and 'inpRays' are part of > # the 'ptsClass' set of classes and are in separate files.) > > # inpRect > # inpScale > > # 27 Jul 2025 / Richard Raubertas > > oldopt <- lsetPkgOpt(reset=TRUE) > #----- > cat(" inpRect ... ") inpRect ... > #== Testing > fspec <- fnSpec(inpnames=c("X1", "X2", "X3")) > # Directly create a hyper-rectangle: > rect1 <- inpRect(rbind(c(0, -1.1, -4), c(2, -0.5, 1000)), spec=fspec) > > # Accessor functions. > corners <- rectCorners(rect1) > stopifnot(dim(corners) == c(2, inpDim(fspec)), + inpNames(rect1) == inpNames(fspec), + colnames(corners) == inpNames(fspec)) > refpt <- rectRefPt(rect1) > stopifnot(dim(refpt) == c(1, inpDim(fspec)), + colnames(refpt) == inpNames(fspec)) > size <- rectSize(rect1) > stopifnot(length(size) == inpDim(fspec), + names(size) == inpNames(fspec), + isTRUE(all.equal(size, corners[2, ] - corners[1, ]))) > > tools::assertError( # hyper-rectangle must be finite + inpRect(rbind(c(0, -1.1, -4), c(2, -0.5, Inf)), spec=fspec), + verbose=warn.tools_verbose) > # Hyper-rectangle must have non-zero d-dimensional volume. > tools::assertError( + inpRect(rbind(c(0, -1.1, 4), c(2, -0.5, 4 + 1e-8)), spec=fspec), + verbose=warn.tools_verbose) > > # Non-default reference point: > rect1 <- inpRect(rbind(c(0, 1, 1), c(0.5, 3, 6)), refpt=c(0, 1, 4), + spec=fspec) > # Update reference point: > rect2 <- update(rect1, refpt=c(0.2, 3, 3)) > stopifnot(isTRUE(all.equal(as.vector(rectRefPt(rect2)), c(0.2, 3, 3)))) > > tools::assertError( # 'refpt' not inside the rectangle + rect2 <- update(rect1, refpt=c(0.2, 0, 3)), + verbose=warn.tools_verbose) > > # If refpt is only slightly outside the rectangle, the rectangle will > # be adjusted to accomodate it. > fspec <- fnSpec(inpnames=c("X1", "X2", "X3"), inptol=c(0, 1e-3)) # defines "slightly" > rect1 <- inpRect(rbind(c(0, 1, 1), c(0.5, 3, 6)), refpt=c(0, 1, 4), + spec=fspec) > rect1b <- update(rect1, refpt=c(0.2, 3 + 5e-4, 3)) > # Strictly includes the new refpt: > stopifnot(inRect(c(0.2, 3 + 5e-4, 3), rect=rect1b, bdryOnly=TRUE, + inptol=inpTol(fspec)/1000), + all(rectCorners(rect1b)[2, ] >= c(0.2, 3 + 5e-4, 3))) > # Not true for original rect: > stopifnot(!all(rectCorners(rect1)[2, ] >= c(0.2, 3 + 5e-4, 3))) > > cat("OK\n") OK > > #----- > cat(" inpScale ... ") inpScale ... > #== Testing > dat <- data.matrix(iris[, 1:4]) > fspec <- fnSpec(inpnames=colnames(dat)) > covmat <- cov(dat) > sds <- sqrt(diag(covmat)) > cormat <- cor(dat) > > # Check equivalence of alternate ways of specifying V. > scl1 <- inpScale(sds, spec=fspec) # w = 1/sds > stopifnot(isTRUE(all.equal(scl1$w, 1/sds)), scl1$V_is_diag, + isTRUE(all.equal(diag(scl1$W), unname(1/sds))), + isTRUE(all.equal(diag(scl1$Winv), unname(sds)))) > scl1a <- inpScale(diag(sds^2, nrow=4), spec=fspec) > scl2 <- inpScale(covmat, spec=fnSpec(scl1)) # 'spec' copied from 'scl1' > scl2a <- inpScale(covmat, spec=fspec) # consistency check > stopifnot(isTRUE(all.equal(scl1$W, scl1a$W)), + identical(scl2, scl2a)) > > # Identity checks. > stopifnot(isTRUE(all.equal(scl2$Winv %*% scl2$W, diag(nrow=ncol(covmat)))), + isTRUE(all.equal(t(scl2$W) %*% scl2$W, unname(solve(covmat)))), + isTRUE(all.equal(scl2$Winv %*% t(scl2$Winv), unname(covmat)))) > > # Mahalanobis distance and Euclidean distance (from origin in this example). > Vinv <- solve(covmat) > mdist <- apply(dat, 1, function(y) { sqrt( t(y) %*% Vinv %*% y) }) > scl2 <- inpScale(covmat, spec=fspec) > tdat <- dat %*% t(scl2$W) # t() since we are transforming _row_ vectors > edist <- sqrt(rowSums(tdat^2)) # Euclidean dist using transformed vectors > stopifnot(isTRUE(all.equal(mdist, edist))) > > # Coordinate system rescaling/rotation. > chk1 <- sweep(dat, 2, scl1$w, "*") # re-scale each dimension to unit SD > stopifnot(isTRUE(all.equal(cov(chk1), cormat))) > chk2 <- dat %*% t(scl2$W) # re-scale and orthogonalize (using t() since > # we are transforming _row_ vectors) > stopifnot(isTRUE(all.equal(cov(chk2), diag(4), check.attributes=FALSE))) > # Reverse the transformation: > z_to_y <- scl2$Winv > chk3 <- chk2 %*% t(z_to_y) # t() since we are transforming _row_ vectors > dimnames(chk3) <- dimnames(dat) > stopifnot(isTRUE(all.equal(chk3, dat))) > > # Empty/null 'inpScale' ('w', a vector of 1's, 'V' the identity matrix): > scl <- inpScale(spec=fspec) > stopifnot(identical(unname(scl$w), rep(1, inpDim(fspec))), + names(scl$w) == inpNames(fspec), + identical(unname(scl$V), diag(nrow=inpDim(fspec))), + identical(dimnames(scl$V), list(inpNames(fspec), inpNames(fspec))), + scl$V_is_diag) > > #----- 1D example. > fspec <- fnSpec(inpnames="X1") > scl <- inpScale(matrix(5, nrow=1, ncol=1), spec=fspec) > stopifnot(identical(fnSpec(scl), fspec)) > > cat("OK\n") OK > > #----- > lsetPkgOpt(oldopt) > > cat('inpRays.r', sep='\n') inpRays.r > > # 15 Feb 2025 / Richard Raubertas > # 25 Feb 2025 : Include tests of new 'update' method, including preservation > # of ray directions. 'rayOrigin' now includes colnames. > > oldopt <- lsetPkgOpt(reset=TRUE) > #----- > cat(" inpRays ... ") inpRays ... > #== Testing > fspec <- fnSpec(inpnames=c("X1", "X2", "X3")) > d <- inpDim(fspec) > > #-- Rays defined by position 1 pts > pts <- rbind(c(0, 1, 3), c(2, 2, 4), c(5, -1, 2)) > orig <- c(0.2, 1.1, 5.5) > rays <- inpRays(pts, origin=orig, spec=fspec) > # Or by specifying their direction vectors and origin: > vecs <- sweep(pts, 2, orig, "-") # equivalent to 'rayVec(rays)' except dimnames > stopifnot(isTRUE(all.equal(vecs, unname(rayVec(rays))))) > rays2 <- inpRays(vecs=vecs, origin=orig, spec=fspec) > stopifnot(identical(rays, rays2)) > > #-- Ray names: > rownames(pts) <- paste0("ray", seq_len(nrow(pts))) > rays1 <- inpRays(pts, spec=fspec) > stopifnot(names(rays1) == rownames(pts)) > vecs <- sweep(pts, 2, orig, "-") > rays2 <- inpRays(vecs=pts, spec=fspec) > stopifnot(identical(names(rays2), names(rays1))) > names(rays1) <- LETTERS[1:3] > stopifnot(names(rays1) == LETTERS[1:3]) # changed > > #-- Subscripting: > stopifnot(identical(rays[], rays)) > stopifnot(length(rays[2:3]) == 2) > stopifnot(length(rays[0]) == 0) > > #-- Accessors and methods. > stopifnot(dim(rayOrigin(rays)) == c(1, d), + dim(rayPosn1(rays)) == c(length(rays), d), + dim(rayVec(rays)) == c(length(rays), d)) > rayrect <- rayRect(rays) > if (!is.null(rayrect)) { + stopifnot(inherits(rayrect, "inpRect"), + identical(rectRefPt(rayrect), rayOrigin(rays)), + all(inRect(rayPosn1(rays), rect=rayrect, bdryOnly=TRUE))) + } > # Can update the origin or lengths of ray vectors with the (internal) > # 'scale' method; see the examples for 'scale.inpRays'. > > rays2 <- inpRays(rbind(c(1, 1, 1), c(1, 2, 3)), origin=c(0, 0, 0), spec=fspec) > # 'rayVecs' and 'rayPosn1' are equivalent when the ray origin is 0: > stopifnot(isTRUE(all.equal(rayVec(rays2), rayPosn1(rays2)))) > > #-- Empty 'inpRays': > rays0 <- inpRays(spec=fspec) # no rays, 0 origin > stopifnot(inherits(rays0, "inpRays"), length(rays0) == 0) > rays0 <- inpRays(spec=fspec, origin=orig) # no rays > stopifnot(length(rays0) == 0, rayOrigin(rays0) == orig) > stopifnot(identical(rays0, rays[0])) > > #-- Combining objects: > stopifnot(identical(rays, c(rays[0], NULL, rays[1:2], rays[3]))) > rays2 <- inpRays(vecs=rayVec(rays), origin=rayOrigin(rays) + 2, spec=fspec) > tools::assertError(rays12 <- c(rays, rays2), # different origins + verbose=warn.tools_verbose) > > stopifnot(identical(rays, inpRays(rays))) > > #-- Ray position 1 points can be (re)defined by each ray's intersection > # with the boundary of a hyper-rectangle: > rect1 <- inpRect(rbind(c(0,0,0), c(1,1,1)), spec=fspec) > # (default reference pt is center: c(0.5, 0.5, 0.5)) > vecs1 <- rbind(c(-1,-1,0), c(0,0,1), c(0, -1, 1)) > rays1 <- inpRays(vecs=vecs1, + rect=rect1) # default 'origin' and 'spec' taken from 'rect' > stopifnot(all(inRect(rayPosn1(rays1), rect=rect1, bdryOnly=TRUE)), + identical(rayOrigin(rays1), rectRefPt(rect1)), + identical(rayRect(rays1), rect1)) > > # Explicit 'origin' argument overrides use of reference point from > # 'rect', and updates the latter: > rays1 <- inpRays(vecs=vecs1, origin=c(1, 1, 0.5), # origin is on boundary of rect + rect=rect1, spec=fspec) > stopifnot(rayOrigin(rays1) == c(1, 1, 0.5)) > stopifnot(inRect(rayPosn1(rays1), rect=rect1, bdryOnly=TRUE)) > stopifnot(identical(rectRefPt(rayRect(rays1)), rayOrigin(rays1))) > > # An existing 'inpRays' object can be updated with a new origin and/or > # hyper-rectangle (which must contain the existing or new ray origin). > # Position 1 points are on the boundary of the new rectangle, the new > # rectangle reference pt is set to the ray origin, which is unchanged. > rect2 <- inpRect(rbind(c(-1, -2, -3), c(1, 2, 3)), spec=fspec) > rays2 <- update(rays1, rect=rect2) > stopifnot(identical(rectCorners(rayRect(rays2)), rectCorners(rect2)), + identical(rayOrigin(rays2), rayOrigin(rays1)), + inRect(rayPosn1(rays2), rect=rect2, bdryOnly=TRUE), + # Ray directions are unchanged: + isTRUE(all.equal(rowSums(rayVec(rays1)*rayVec(rays2)), + sqrt(rowSums(rayVec(rays1)^2) * + rowSums(rayVec(rays2)^2))))) > # Note: internally, ref pt of 'rect' is set to ray origin: > stopifnot(identical(rectRefPt(rayRect(rays2)), rayOrigin(rays2))) > # Change only 'origin' (which should also reset ref pt of its rect). > rays3 <- update(rays1, origin=c(0.2, 0.3, 0.4)) > stopifnot(rectCorners(rayRect(rays3)) == rectCorners(rect1), + as.vector(rayOrigin(rays3)) == c(0.2, 0.3, 0.4), + as.vector(rectRefPt(rayRect(rays3))) == c(0.2, 0.3, 0.4)) > # Change both 'origin' and 'rect'. > rays3 <- update(rays1, origin=c(0, -1, 2), rect=rect2) > stopifnot(identical(rectRefPt(rayRect(rays3)), rayOrigin(rays3)), + inRect(rayPosn1(rays3), rect=rect2, bdryOnly=TRUE), + as.vector(rayOrigin(rays3)) == c(0, -1, 2), + # Ray directions are unchanged: + isTRUE(all.equal(rowSums(rayVec(rays1)*rayVec(rays3)), + sqrt(rowSums(rayVec(rays1)^2) * + rowSums(rayVec(rays3)^2))))) > > # Rays with different rect's can be combined, as long as they have the > # same origin. The combination uses the rect of the first. > rays12 <- c(rays1, rays2) > stopifnot(length(rays12) == length(rays1) + length(rays2), + identical(rayRect(rays12), rayRect(rays1))) > > #-- Degenerate rays are dropped: > vecs2 <- rbind(c(-1,-1,0), c(0,0,0), c(0, -1, 1)) # 2nd vector has 0 length > tools::assertWarning( + rays <- inpRays(vecs=vecs2, origin=c(0.5, 0.5, 0.5), spec=fspec), + verbose=warn.tools_verbose) > stopifnot(length(rays) == 2) # not 3 > # One can also end up with a degenerate ray if a rect is specified and > # its reference point is on the boundary. In that case a ray's only > # intersection with the rect boundary may be the origin: > rect <- inpRect(rbind(c(0,0,0), c(1,1,1)), spec=fspec) > orig <- c(1,1,0.5) # on the boundary of 'rect' > vecs <- rbind(c(-1,-1,0), c(0,1,1), c(0, -1, 1)) > # The second direction vector points away from the interior of 'rect', > # relative to 'orig'. > tools::assertWarning( + rays <- inpRays(vecs=vecs, origin=orig, rect=rect), + # (Second ray's only intersection with 'rect' is the origin.) + verbose=warn.tools_verbose) > stopifnot(length(rays) == 2) > # Can force degenerate rays to be treated as an error: > tools::assertError( + rays <- inpRays(vecs=vecs, origin=orig, rect=rect, warnDegenerate=2), + verbose=warn.tools_verbose) > > #-- Other warnings and errors. > orig2 <- c(2, 2, 6) > tools::assertError( + inpRays(pts, origin=orig2, rect=rect1), + # ray origin does not belong to rectangle + verbose=warn.tools_verbose) > > cat("OK\n") OK > > #----- > lsetPkgOpt(oldopt) > cat('inpRaysMethods.r', sep='\n') inpRaysMethods.r > # Methods specifically for 'inpRays' objects (where corresponding 'inpPts' > # and 'inpPairs' methods in 'ptsClassMethods.r' are not adequate). > # > # scale.inpRays* (16 Jan 2026: works ok, has tests, but not currently > # used, so disabled) > # plot.inpRays > # update.inpRays > > # *= not public > > # 13 Apr 2026 / Richard Raubertas > > oldopt <- lsetPkgOpt(reset=TRUE) > #----- > if (exists("scale.inpRays")) { + cat(" scale.inpRays ... ") + #== Testing + fspec <- fnSpec(inpnames=LETTERS[1:3]) + pRect <- inpRect(rbind(c(2, 0, -2), c(4, 4, 4)), + spec=fspec) # 2 x 4 x 6 box centered at (3,2,1) + pRays <- inpRays(rbind(c(2, 0, -2), + c(4, 0, 1), + c(3, 2, 3), + c(10, 10, 10)), + origin=c(3, 3, 3), spec=fspec) + scRays1 <- scale(pRays, center=TRUE, scale=FALSE) + scRays2 <- scale(pRays, center=c(3, 2, 1), scale=FALSE) + scRays3 <- scale(pRays, center=pRect, scale=FALSE) + stopifnot(rayOrigin(scRays1) == t(rep(0, 3)), + identical(scRays2, scRays3)) + + chk <- rayVec(scale(pRays, center=FALSE, scale=TRUE)) # unit length + stopifnot(isTRUE(all.equal(sqrt(rowSums(chk^2)), rep(1, nrow(chk))))) + + # Translate and scale ray vectors to match reference point and boundary + # of a hyper-rectangle: + scRays4 <- scale(pRays, center=c(3, 2, 1), scale=pRect) + scRays5 <- scale(pRays, center=pRect, scale=pRect) + stopifnot(identical(scRays4, scRays5)) + raypts5 <- rayPosn1(scRays5) + stopifnot(all(inBox(raypts5, corners=rectCorners(pRect), bdryOnly=TRUE, + inptol=inpTol(fspec)))) + scRays6 <- scale(pRays, center=FALSE, scale=pRect) + raypts6 <- rayPosn1(scRays6) + stopifnot(all(inBox(raypts6, corners=rectCorners(pRect), bdryOnly=TRUE, + inptol=inpTol(fspec)))) + + # Scale lengths of ray vectors without changing origin. + scRays7 <- scale(pRays, center=FALSE, scale=2) # *divides* by 2 + len <- sqrt(rowSums(rayVec(pRays)^2)) + len7 <- sqrt(rowSums(rayVec(scRays7)^2)) + stopifnot(isTRUE(all.equal(len7, len/2))) + factors <- c(1, 3, 0.5, -2) # negative scale reverses ray vector direction + scRays8 <- scale(pRays, center=FALSE, scale=factors) # *divides* by 'factors' + len <- sqrt(rowSums(rayVec(pRays)^2)) + len8 <- sqrt(rowSums(rayVec(scRays8)^2)) + stopifnot(isTRUE(all.equal(len8, len/abs(factors)))) + # Not allowed to use scale to create infinite or degenerate ray vectors: + factors[c(1, 3)] <- c(0, Inf) + tools::assertError(scale(pRays, center=FALSE, scale=factors), + verbose=warn.tools_verbose) + + #----- 2D example for rotation/scaling of ray vectors by a matrix. + + # Create ray vectors distributed along the ellipse + # X1^2 + 2*X1*X2 + 2*(X2^2) = 8, + # then translate to have origin at c(2, 1): + fspec <- fnSpec(inpnames=c("X1", "X2")) + X1 <- c(-4:4) + X2 <- c((-X1 - sqrt(16 - X1^2))/2, (-X1 + sqrt(16 - X1^2))/2) + X1 <- c(X1, X1) + pRays <- inpRays(vecs=cbind(X1, X2), origin=c(2, 1), spec=fspec) + # Plot them: + #plot(pRays, eqAxes=TRUE, rayPar=list(col="blue")) + + # The untranslated ellipse can also be expressed as points X with a + # Mahalanobis distance of 8 from the origin, given by the quadratic + # form t(X) %*% V^(-1) %*% X, where symmetric, positive-definite matrix + # V is + V <- matrix(c(2, -1, -1, 1), nrow=2) + stopifnot(isTRUE(all.equal(diag(cbind(X1, X2) %*% solve(V) %*% rbind(X1, X2)), + rep(8, length(pRays))))) # all 8 + + # An 'inpScale' object that transforms this ellipse to a unit circle is: + pScale <- inpScale(8*V, spec=fspec) + # ('8*' because ultimately we want to divide vectors by sqrt(8) to get + # unit length--'mat' acts quadratically not linearly on coordinates.) + # Alternative: inpScale(x=sqrt(c(8, 8)), mat=V, spec=fspec) + scRays1 <- scale(pRays, center=FALSE, scale=pScale) + stopifnot(isTRUE(all.equal(rowSums(rayVec(scRays1)^2), + rep(1, length(pRays))))) # unit length + #plot(scRays1, add=TRUE, rayPar=list(col="red")) + # All have unit length, but are not in the same directions as the original + # vectors (rotated). To keep the same directions, instead scale each + # vector individually by its length: + raylens <- sqrt(rowSums(rayVec(pRays)^2)) + scRays2 <- scale(pRays, center=FALSE, scale=raylens) + stopifnot(isTRUE(all.equal(rowSums(rayVec(scRays2)^2), + rep(1, length(pRays))))) # unit length + #plot(scRays2, add=TRUE, rayPar=list(col="green")) + + cat("OK\n") + } > > #----- > cat(" update.inpRays ... ") update.inpRays ... > fspec <- fnSpec(inpnames=c("X1", "X2")) > set.seed(1) > rays <- randomRays(10, spec=fspec) # ray vectors have unit length, > # origin at (0, 0) > rect2 <- inpRect(rbind(c(0, 0), c(1, 1)), spec=fspec) # unit square > rays2 <- update(rays, origin=c(0.3, 0.4), rect=rect2) > stopifnot(identical(as.vector(rayOrigin(rays2)), c(0.3, 0.4)), + identical(rectCorners(rayRect(rays2)), rectCorners(rect2)), + # Ref pt has been set to new origin: + identical(rectRefPt(rayRect(rays2)), rayOrigin(rays2)), + inRect(rayPosn1(rays2), rect=rect2, bdryOnly=TRUE)) > > rays3 <- update(rays2, origin=c(0.5, 0.6)) > stopifnot(identical(as.vector(rayOrigin(rays3)), c(0.5, 0.6)), + identical(rectCorners(rayRect(rays3)), rectCorners(rect2)), + # Ref pt has been set to new origin: + identical(rectRefPt(rayRect(rays3)), rayOrigin(rays3)), + inRect(rayPosn1(rays3), rect=rect2, bdryOnly=TRUE)) > cat("OK\n") OK > > #----- > lsetPkgOpt(oldopt) > cat('lsetPkgOpt.r', sep='\n') lsetPkgOpt.r > # User-level and internal package options. > # > # lsetPkgOpt > # internalOpt > # < create environment 'lvlSets_env'> > # .onLoad > > # 09 Feb 2025 / Richard Raubertas > # 17 Feb 2026 : Change 'tol.resp', 'tol.ptCoord' to 'resptol', 'inptol', for > # consistency with other functions. > > oldopt <- lsetPkgOpt(reset=TRUE) > #----- > cat(" lsetPkgOpt ... ") lsetPkgOpt ... > #== Testing > oldopts0 <- lsetPkgOpt(reset=TRUE) # don't assume current opts are default > # Default 'factory-fresh' options: > defopt <- lsetPkgOpt() > #str(defopt) > > tol1 <- lsetPkgOpt("resptol") # plain vector > stopifnot(is.vector(tol1), length(tol1) == 2) > oldtol1 <- lsetPkgOpt(resptol=1e-4) # list of length 1 > stopifnot(identical(tol1, oldtol1[[1]])) > > opts <- lsetPkgOpt(c("warn.invalidProfile", "warn.showBadRespPts")) > stopifnot(is.list(opts), length(opts) == 2) > oldopts <- lsetPkgOpt(warn.invalidProfile=0, warn.showBadRespPts=10) > stopifnot(identical(opts, oldopts)) > lsetPkgOpt(reset=TRUE) > stopifnot(identical(lsetPkgOpt(), defopt)) > oldopts2 <- lsetPkgOpt(list(warn.invalidProfile=0, warn.showBadRespPts=10), + foo="bar") # 'foo' ignored > stopifnot(identical(oldopts, oldopts2)) > > # Errors from non-existant options. > tools::assertError(lsetPkgOpt(c("inptol", "foo")), + verbose=warn.tools_verbose) > tools::assertError(lsetPkgOpt(inptol=c(2e-5, 3e-8), foo="bar"), + verbose=warn.tools_verbose) > # Error from setting an invalid value. > tools::assertError(lsetPkgOpt(inptol=c(2e-5, 3e-8), + resptol=c(1e-1)), verbose=warn.tools_verbose) > > # Return to options where we started. > lsetPkgOpt(oldopts0) > stopifnot(identical(lsetPkgOpt(), oldopts0)) > > cat("OK\n") OK > > #----- > cat(" internalOpt ... ") internalOpt ... > #== Testing > oldopts0 <- internalOpt(reset=TRUE) # don't assume current opts are default > # Default 'factory-fresh' options: > defopt <- internalOpt() > #str(defopt) > stopifnot(is.list(oldopts0), names(oldopts0) == names(defopt)) > > oldval <- internalOpt(tol.line_position=c(1e-6, 1.5e-8)) # list w/ one comp > stopifnot(is.list(oldval), length(oldval) == 1) > chk <- internalOpt("tol.line_position") # new value (vector) > stopifnot(!is.list(chk), identical(chk, c(1e-6, 1.5e-8))) > internalOpt(reset=TRUE) > stopifnot(identical(defopt, internalOpt())) > > tools::assertError(internalOpt("foo"), verbose=warn.tools_verbose) # non-existant option > tools::assertError(internalOpt(foo=42), verbose=warn.tools_verbose) # non-existant option > tools::assertError(internalOpt(tol.line_position=c(1e-2, 1e-8)), + verbose=warn.tools_verbose) # invalid value > tools::assertError(internalOpt(tol.line_position=c(1e-2, 1e-8)), + verbose=warn.tools_verbose) # invalid value > > # Restore starting opts. > internalOpt(oldopts0) > > cat("OK\n") OK > > #----- > lsetPkgOpt(oldopt) > cat('lsetSegs.r', sep='\n') lsetSegs.r > # lsetSegs (S3 generic) > # lsetSegs.default* (get '$segs' component from a list) > # lsetSegs.lsetSegs* (no-op, but needed to avoid wrong method dispatch) > # lsetSegs.respProfiles* > # lsetSegs.respPairs* (create from 'respPairs' by adding 'inpLines') > # lsetSegs.fnSpec* (create an empty object) > > # *= not public > > oldopt <- lsetPkgOpt(reset=TRUE) > > #----- > cat(" lsetSegs.default ... ") lsetSegs.default ... > Ex <- dbl_ellipse_X1bnd_2dEx # see '?dbl_ellipse_2dEx' > fobj <- Ex$fobj > segs <- Ex$segs > chk <- lsetSegsCheck(segs, fnobj=fobj, posns=(1:4)/5) # list > stopifnot(identical(lsetSegs(chk), Ex$segs_checked)) > > cat("OK\n") OK > > #----- > cat(" lsetSegs.lsetSegs ... ") lsetSegs.lsetSegs ... > # Should be a no-op. > Ex <- dbl_ellipse_X1bnd_2dEx # see '?dbl_ellipse_2dEx' > segs <- Ex$segs > stopifnot(identical(segs, lsetSegs(segs))) > > cat("OK\n") OK > > #----- > cat(" lsetSegs.fnSpec ... ") lsetSegs.fnSpec ... > #== Testing > Ex <- dbl_ellipse_2dEx > respF <- Ex$fobj > fspec <- fnSpec(respF) > thresh <- Ex$thresh > lsegs <- lsetSegs(fspec, lsetthresh=thresh, checkArg=TRUE) # empty > stopifnot(inherits(lsegs, "lsetSegs"), length(lsegs) == 0) > > # Non-empty 'lines', but no segments. > lines <- axisRays(fspec) > lsegs <- lsetSegs(fspec, lsetthresh=thresh, lines=lines, checkArg=TRUE) # empty > stopifnot(inherits(lsegs, "lsetSegs"), length(lsegs) == 0) > > cat("OK\n") OK > > #----- > lsetPkgOpt(oldopt) > cat('lsetSegsCheck.r', sep='\n') lsetSegsCheck.r > #----- > cat(" lsetSegsCheck ... ") lsetSegsCheck ... > #== Testing > Ex <- dbl_ellipse_X1bnd_2dEx > segs <- Ex$segs # 67 segs on 50 rays, 8 invalid segs > chk <- lsetSegsCheck(segs, fnobj=Ex$fobj, posns=(1:4)/5) > stopifnot(inherits(chk$segs, "lsetSegs"), + identical(as.vector(table(chk$validIn)), c(8L, 59L)), + identical(as.vector(table(chk$validOut)), c(75L))) > > # Should be able to handle an empty 'segs'. > chk <- lsetSegsCheck(segs[0], fnobj=Ex$fobj, posns=(1:4)/5) > stopifnot(inherits(chk$segs, "lsetSegs"), length(chk$segs) == 0, + identical(as.vector(table(chk$validIn)), integer(0)), + identical(as.vector(table(chk$validOut)), integer(0))) > > cat("OK\n") OK > > #----- > cat('lsetSegsMethods.r', sep='\n') lsetSegsMethods.r > # Methods for 'lsetSegs' objects. (Most standard generics use the '.inpPts' > # or '.inpPairs' methods for these objects.) > # > # plot.lsetSegs > # pairs.lsetSegs > # summary.lsetSegs > # print.summary.lsetSegs* > # plotSegsList (not strictly a method, but designed for 'lsetSegs' objects) > > # *= not public > > oldopt <- lsetPkgOpt(reset=TRUE) > #----- > cat(" summary.lsetSegs ... ") summary.lsetSegs ... > #== Testing > # Test that function works even with empty input. > Ex <- dbl_ellipse_2dEx > respF <- Ex$fobj > fspec <- fnSpec(respF) > thresh <- Ex$thresh > lsegs <- lsetSegs(fspec, lsetthresh=thresh, checkArg=TRUE) # empty > summ <- summary(lsegs) > stopifnot(summ$nlines == 0, dim(summ$bBox) == c(2, 2), all(is.na(summ$bBox)), + length(summ$respRange) == 2, all(is.na(summ$respRange)), + is.data.frame(summ$lineInfo), nrow(summ$lineInfo) == 0) > > # Non-empty 'lines', but no segments. > lines <- axisRays(fspec) > lsegs <- lsetSegs(fspec, lsetthresh=thresh, lines=lines, checkArg=TRUE) # empty > summ <- summary(lsegs) > stopifnot(summ$nlines == length(lines), dim(summ$bBox) == c(2, 2), + all(is.na(summ$bBox)), + length(summ$respRange) == 2, all(is.na(summ$respRange)), + is.data.frame(summ$lineInfo), nrow(summ$lineInfo) == length(lines)) > > cat("OK\n") OK > > #----- > lsetPkgOpt(oldopt) > cat('respAccessors.r', sep='\n') respAccessors.r > # Functions to extract response-related components or properties of objects > # in the 'ptsClass' set of classes. > > # lsetThresh > # respInfo (generic and 'respPts'*, 'respPairs'*, 'lsetSegs', 'default' > # methods) > > # *= not public > > # 24 Feb 2025 / Richard Raubertas > # 12 Feb 2026 : Update 'lsetThresh' tests. > > oldopt <- lsetPkgOpt(reset=TRUE) > #----- > cat(" lsetThresh ... ") lsetThresh ... > stopifnot(is.na(lsetThresh(NULL))) > > Ex <- dbl_ellipse_2dEx > segs <- Ex$segs > stopifnot(identical(lsetThresh(segs), Ex$thresh)) > > cat("OK\n") OK > > #----- > cat(" respInfo ... ") respInfo ... > respF <- function(x, inclGrad=FALSE) { + resp <- x[, "X1"] + x[, "X2"]^2 + x[, "X3"]^3 + if (inclGrad) { + grad <- cbind(1, 2*x[, "X2"], 3*x[, "X3"]^2) + structure(resp, gradient=grad) + } else resp + } > fnobj <- fnObj(c("X1", "X2", "X3"), respfn=respF, hasgrad=TRUE) > > # NA point coordinates are not allowed. > tools::assertError(respInfo(c(NA, 1, 3), fnobj=fnobj, inclGrad=TRUE), + verbose=warn.tools_verbose) > > # Evaluate the response function at a set of points. > pts <- rbind(c(0, 1, 3), c(2, 2, 4), c(5, -1, 2)) # three input space points > rpts <- respInfo(pts, fnobj=fnobj, lsetthresh=28) > stopifnot(is.data.frame(rpts), nrow(rpts) == nrow(pts), + is.null(attr(rpts, "gradient"))) > rpts <- respInfo(pts, fnobj=fnobj, lsetthresh=28, inclGrad=TRUE) > grad <- attr(rpts, "gradient") > stopifnot(is.data.frame(rpts), nrow(rpts) == nrow(pts), + is.matrix(grad), dim(grad) == dim(pts), + identical(dimnames(grad), list(NULL, inpNames(fnobj)))) > > # If 'fnobj' supports gradients but they were not actually calculated, > # return all NA's. > rpts <- respPts(pts, fnobj=fnobj, lsetthresh=28, inclGrad=FALSE) > rinfo <- respInfo(rpts, inclGrad=TRUE) > grad <- attr(rinfo, "gradient") > stopifnot(is.data.frame(rinfo), nrow(rinfo) == nrow(pts), + is.matrix(grad), dim(grad) == dim(pts), + identical(dimnames(grad), list(NULL, inpNames(fnobj))), + all(is.na(grad))) > > # Empty input. > resp <- respInfo(pts[0, , drop=FALSE], fnobj=fnobj, inclGrad=TRUE) > stopifnot(is.data.frame(resp), nrow(resp) == 0, + is.matrix(attr(resp, "gradient")), + nrow(attr(resp, "gradient")) == 0) > > cat("OK\n") OK > > #----- > lsetPkgOpt(oldopt) > cat('respProfiles.r', sep='\n') respProfiles.r > # respProfiles (S3 generic) > # (no default method) > # respProfiles.fnObj (calculate from scratch) > # respProfiles.lsetSegs* (extract from an 'lsetSegs' object) > > # *= not public > > # 15 Feb 2025 / Richard Raubertas > # 24 Feb 2025 : 'fnObj' no longer has 'spec' arg. > # 16 Jul 2025 : Take advantage of, and test, 'update' method for 'fnObj'. > # 13 Feb 2026 : 3D and 5D examples moved to 'tests_extra'. > > oldopt <- lsetPkgOpt(reset=TRUE) > #----- > cat(" respProfiles.fnObj ... \n") respProfiles.fnObj ... > #== Testing > #----- 1D example. > cat(" 1D example ... ") 1D example ... > respf <- function(x, ...) { x*cos(x) } > fobj <- fnObj(c("X1"), respfn=respf) # default tolerances > ray1 <- inpRays(c(1), spec=fobj) # unit vector > rslt <- respProfiles(fobj, rays=ray1, posns=seq(-3, 10), lsetthresh=0.2) > > # Methods for 'respProfiles' objects: > summ <- summary(rslt) > stopifnot(length(rslt) == 1, is.list(summ), identical(inpLines(rslt), ray1), + identical(summ$lsetThresh, lsetThresh(rslt)), + is.matrix(summ$ptCounts), dim(summ$ptCounts) == c(summ$nlines, 5), + colSums(summ$ptCounts) == c(19, 19, 0, 11, 5), + is.matrix(summ$respRange), dim(summ$respRange) == c(summ$nlines, 2)) > > # Extract and summarize level set segments. > lssegs <- lsetSegs(rslt) # '.respProfiles' method > summ <- summary(lssegs) > stopifnot(length(lssegs) == 3, is.list(summ), + identical(inpLines(lssegs), ray1), + identical(summ$lsetThresh, lsetThresh(lssegs)), + is.matrix(summ$bBox), dim(summ$bBox) == c(2, inpDim(lssegs)), + is.data.frame(summ$lineInfo), nrow(summ$lineInfo) == length(ray1), + colSums(summ$lineInfo[, 1:4]) == c(length(lssegs), 0, 5, 1), + summ$nlines == nrow(summ$lineInfo), length(summ$respRange) == 2) > > rpts <- respPts(rslt) > summ <- summary(rpts) > stopifnot(length(rpts) == 19, is.list(summ), + identical(summ$lsetThresh, lsetThresh(rpts)), + is.matrix(summ$ptCounts), dim(summ$ptCounts) == c(1, 3), + colSums(summ$ptCounts) == c(19, 19, 11), + length(summ$respRange) == 2) > > # Edge case > summ <- summary(rpts[0]) > stopifnot(is.list(summ), + identical(summ$lsetThresh, lsetThresh(rpts)), + is.matrix(summ$ptCounts), dim(summ$ptCounts) == c(1, 3), + colSums(summ$ptCounts) == c(0, 0, 0), + length(summ$respRange) == 2) > > # Implicit definition of feasible region (response function returns NA > # to indicate infeasible points). Note that feasible region is closed. > respf2 <- function(x, ...) { x <- as.vector(x) + ifelse(x <= 1 | x >= 2, x*cos(x), NA) } > fobj2 <- fnObj(inpNames(fobj), respfn=respf2, warnNAresp=FALSE) > rslt2 <- respProfiles(fobj2, rays=ray1, lsetthresh=0.2, + posns=c(seq(-3, 9), 1.3)) # need a test position in > # the infeasible region to detect it > rinfo2 <- respInfo(rslt2) > stopifnot(nrow(rinfo2) == 18, + rinfo2$feas == c(rep(TRUE, 7), FALSE, rep(TRUE, 10)), + identical(rinfo2$feasbdry, c(NA, rep(FALSE, 5), TRUE, FALSE, TRUE, + rep(FALSE, 8), NA)), + rinfo2$lset == c(rep(TRUE, 3), rep(FALSE, 2), rep(TRUE, 2), + rep(FALSE, 4), rep(TRUE, 5), rep(FALSE, 2)), + identical(rinfo2$lsetbdry, c(NA, FALSE, TRUE, FALSE, FALSE, TRUE, + rep(FALSE, 5), TRUE, rep(FALSE, 3), + TRUE, FALSE, FALSE))) > cat("OK\n") OK > > #----- Multimodal example in 2D (surface with two elliptical contours). > cat(" 2D multimodal example ... ") 2D multimodal example ... > respf <- function(x) { + c1 <- x[, 1]^2 + 2*x[, 1]*x[, 2] + 2*(x[, 2]^2) + # Same ellipse, but shifted to be centered at (8, 8): + c2 <- (x[, 1]-8)^2 + 2*(x[, 1]-8)*(x[, 2]-8) + 2*((x[, 2]-8)^2) + resp <- -1 * pmin(c1, c2) + resp + } > # Feasible region: X1 >= 0. > feasf <- function(x) { x <- num_matrix(x); (x[, 1] >= 0) } > fobj <- fnObj(c("X1", "X2"), respfn=respf, feasfn=feasf) > > # Calculate ray profiles in four standard directions. > rays1 <- axisRays(fobj, degree=2) > rslt1 <- respProfiles(fobj, rays=rays1, posns=seq(-8, 15), lsetthresh=-40) > summ <- summary(rslt1) > stopifnot(length(rslt1) == length(rays1), identical(inpLines(rslt1), rays1), + is.list(summ), identical(summ$lsetThresh, lsetThresh(rslt1)), + is.matrix(summ$ptCounts), dim(summ$ptCounts) == c(summ$nlines, 5), + colSums(summ$ptCounts) == c(101, 77, 3, 48, 8), + is.matrix(summ$respRange), dim(summ$respRange) == c(summ$nlines, 2)) > > # '[' and 'c' methods test. > stopifnot(identical(rslt1, c(rslt1[1], rslt1[0], NULL, rslt1[2:4], NULL))) > > # Extract level set segments. > lssegs <- lsetSegs(rslt1) # .respProfiles method > summ <- summary(lssegs) > stopifnot(length(lssegs) == 6, is.list(summ), + identical(inpLines(lssegs), rays1), + identical(summ$lsetThresh, lsetThresh(lssegs)), + is.matrix(summ$bBox), dim(summ$bBox) == c(2, inpDim(lssegs)), + is.data.frame(summ$lineInfo), nrow(summ$lineInfo) == length(rays1), + colSums(summ$lineInfo[, 1:4]) == c(length(lssegs), 3, 8, 1), + summ$nlines == nrow(summ$lineInfo), length(summ$respRange) == 2) > # Test that 'lsetbdry', 'feasbdry', 'lsetcens' are consistent. (Regression test) > rinfo <- respInfo(lssegs) > stopifnot(is.data.frame(rinfo), nrow(rinfo) == 2*length(lssegs), + !anyNA(rinfo[, "feas"]), !anyNA(rinfo[, "lset"]), + !any(rinfo[!rinfo$feas, "feasbdry"]), + all(rinfo$feas[naF(rinfo$feasbdry)]), + !any(rinfo[!rinfo$lset, "lsetbdry"]), + all(rinfo$lset[naF(rinfo$lsetbdry)]), + !anyNA(rinfo[, "lsetcens"]), + all(rinfo$lsetcens == !(naF(rinfo$feasbdry) | naF(rinfo$lsetbdry)))) > > # Just feas bdry search, no level set. > rslt2 <- respProfiles(fobj, rays=rays1, posns=seq(-8, 15), feasbdrySrch=TRUE) > summ <- summary(rslt2) # same 'nfeasbdry' as for 'rslt' > stopifnot(length(rslt2) == length(rays1), identical(inpLines(rslt2), rays1), + is.list(summ), identical(summ$lsetThresh, lsetThresh(rslt2)), + is.matrix(summ$ptCounts), dim(summ$ptCounts) == c(summ$nlines, 5), + identical(unname(colSums(summ$ptCounts)), c(96, 72, 3, NA, NA)), + is.matrix(summ$respRange), dim(summ$respRange) == c(summ$nlines, 2)) > > # Profile the lines that contain the rays. (Mostly equivalent.) > lines1 <- inpLines(rays1) > rslt1a <- respProfiles(fobj, lines=lines1, posns=seq(-8, 15), lsetthresh=-40) > stopifnot(length(rslt1a) == length(rslt1), + identical(summary(rslt1)$ptCounts, summary(rslt1a)$ptCounts)) > > # Edge cases: no rays; no posns. > rslt3 <- respProfiles(fobj, rays=rays1[0], posns=seq(-8, 15), + lsetthresh=-40) # no rays, no points > summ <- summary(rslt3) > stopifnot(inherits(rslt3, "respProfiles"), length(rslt3) == 0, + nrow(summ$ptCounts) == 0) > rslt4 <- respProfiles(fobj, rays=rays1, posns=numeric(0), + lsetthresh=-40) # 4 rays, no points > summ <- summary(rslt4) > stopifnot(inherits(rslt4, "respProfiles"), length(rslt4) == 4, + nrow(summ$ptCounts) == 4, all(summ$ptCounts == 0)) > stopifnot(identical(rslt4, c(rslt4[1], NULL, rslt4[2:4], rslt4[0]))) > if (requireNamespace("ggplot2", quietly=TRUE)) { + tools::assertError(print(plot(rslt3)), verbose=warn.tools_verbose) # nothing to plot + tools::assertError(print(plot(rslt4)), verbose=warn.tools_verbose) # nothing to plot + } > cat("OK\n") OK > > #----- > lsetPkgOpt(oldopt) > cat('respPts.r', sep='\n') respPts.r > # respPts* (S3 generic) > # respPts.default* (handles matrices and 'inpPts' objects) > # respPts.respPts* (extract or update an existing object. 'strictInfoCols' arg) > # respPts.data.frame* (assemble 'respPts' manually. not public.) > # respPts.fnSpec* (create an empty object) > # respPts.respPairs* (includes 'which' arg. Also handles 'lsetSegs') > > # *= not public > > # 12 Feb 2025 / Richard Raubertas > # 24 Feb 2025 : 'fnObj' no longer has 'spec' arg. > > oldopt <- lsetPkgOpt(reset=TRUE) > > #----- > cat(" respPts ... ") respPts ... > #== Testing > respF <- function(x, inclGrad=FALSE) { + resp <- x[, "X1"] + x[, "X2"]^2 + x[, "X3"]^3 + if (inclGrad) { + grad <- cbind(1, 2*x[, "X2"], 3*x[, "X3"]^2) + structure(resp, gradient=grad) + } else resp + } > fnobj <- fnObj(c("X1", "X2", "X3"), respfn=respF, hasgrad=TRUE) > coord <- rbind(c(0, 1, 3), c(2, 2, 4), c(5, -1, 2)) # three input space points > # 'respPts' object created directly from a response function: > rpts <- respPts(coord, fnobj=fnobj, lsetthresh=28, inclGrad=TRUE) > > # Construct the same 'respPts' from an 'inpPts' object. > ipts <- inpPts(coord, spec=fnobj) > rpts1b <- respPts(ipts, fnobj=fnobj, lsetthresh=28, inclGrad=TRUE) > stopifnot(identical(rpts1b, rpts)) > > # Construct the same 'respPts' object starting with a data frame of > # response function values: > info <- respInfo(rpts, inclGrad=TRUE) # data frame w/ attribute > stopifnot(is.matrix(attr(info, "gradient")), + dim(attr(info, "gradient")) == dim(coord)) > rpts2 <- respPts(info, spec=fnSpec(fnobj), ptcoord=coord, lsetthresh=28, + respgrad=attr(info, "gradient")) > stopifnot(identical(rpts, rpts2)) > rpts2 <- respPts(info, spec=fnSpec(fnobj), ptcoord=coord, lsetthresh=28) > stopifnot(identical(rpts, rpts2)) # tests default for 'respgrad' > # Test code branch that skips most arg checking: > rpts2b <- respPts(info, spec=fnSpec(fnobj), ptcoord=coord, lsetthresh=28, + respgrad=attr(info, "gradient"), checkArg=FALSE) > stopifnot(identical(rpts2, rpts2b)) > > # If one only wants the response values, the data frame can > # be created directly, rather than going through 'respPts': > respinfo <- respInfo(coord, fnobj=fnobj, lsetthresh=28) > stopifnot(identical(respinfo, respInfo(rpts))) > > # Creating an empty 'respPts'. > rpts0 <- respPts(matrix(0, nrow=0, ncol=inpDim(fnobj)), fnobj=fnobj) > stopifnot(inherits(rpts0, "respPts"), length(rpts0) == 0) > rpts0b <- respPts(fnSpec(fnobj)) > stopifnot(identical(rpts0, rpts0b)) > > # Indexing. > stopifnot(length(rpts[]) == length(rpts), length(rpts[2]) == 1, + length(rpts[0]) == 0) > > # Combine multiple 'respPts' objects. > stopifnot(identical(rpts, c(rpts[0], rpts[1], NULL, rpts[2:3]))) > > # Accessor functions. > stopifnot(inherits(fnSpec(rpts), "fnSpec"), + identical(unname(ptCoord(rpts)), coord)) > > # Updating 'lsetThresh', 'respTol'. > rpts2 <- respPts(rpts, lsetthresh=28 + 1e-6) > stopifnot(identical(lsetThresh(rpts2), c(28 + 1e-6)), + respInfo(rpts2)$lset == c(TRUE, TRUE, FALSE)) > # response value 28 is within tolerance of lsetthresh > rpts2b <- respPts(rpts2, resptol=c(1e-8, 1e-8)) > stopifnot(respInfo(rpts2b)$lset == c(FALSE, TRUE, FALSE)) > # with tighter tolerance, response value 28 is no longer > # considered to be in the level set > > # Can't combine objects with different level set definitions: > tools::assertError(c(rpts, rpts2), verbose=warn.tools_verbose) > tools::assertError(c(rpts2, rpts2b), verbose=warn.tools_verbose) > # Can reverse updates: > rpts3 <- respPts(rpts2, lsetthresh=28) # reverse update > stopifnot(identical(rpts, rpts3)) > rpts3b <- respPts(rpts2b, resptol=respTol(rpts)) > stopifnot(identical(rpts2, rpts3b)) > > stopifnot(class(rpts) == c("respPts", "inpPts"), + class(inpPts(rpts)) == "inpPts") > stopifnot(identical(rpts, respPts(rpts))) > > cat("OK\n") OK > > #----- > lsetPkgOpt(oldopt) > > proc.time() user system elapsed 18.81 0.32 19.12