Nothing
##================================================
## Copyright (C) 2003-2018 Fisheries & Oceans Canada
## Nanaimo, British Columbia
## This file is part of PBS Mapping.
##================================================
##------------------------------------------------
## Dot (hidden) functions kept in one file for
## convenience (RH 2018-08-24)
##------------------------------------------------
## .addAxis-----------------------------2011-12-20
## 'xlim'/'ylim': vectors of length 2
## 'tckLab'/'tck'/'tckMinor': vectors of length 1 or 2
## Returns: NULL (invisible)
## ---------------------------------------------NB
.addAxis <- function(xlim, ylim, tckLab, tck, tckMinor, ...)
{
# force 'tck*' to length 2
tckLab <- rep(tckLab, length.out = 2);
tck <- rep(tck, length.out = 2);
tckMinor <- rep(tckMinor, length.out = 2);
# reduce 'cex' and 'mex' to appropriate sizes for axis labels; these settings
# should remain after the function termimates; .addLabels() will use the same
# units for measuring lines
# changing 'cex' causes plot(...) to change 'mai', so let's reset 'mai'
# after setting 'cex'
mai <- par()$mai;
par(cex = par()$cex * 0.8); # decrease font size
par(mai = mai);
par(mex = par()$cex); # decrease spacing to match font size
# 1 is the horizontal axis, 2 is the vertical axis
lim <- list(xlim, ylim);
rotate <- list(0, 90)
for (i in 1:2) {
if (((i == 1) && (par()$xaxt != "n")) ||
((i == 2) && (par()$yaxt != "n"))) {
# create both major and minor ticks
ticks <- pretty(lim[[i]]);
ticksMinor <- pretty(c(0, diff(ticks)[1]));
ticksMinor <- (sort(rep(ticks, length(ticksMinor))) +
rep(ticksMinor, length(ticks)));
# filter them for anything on the extents
ticks <- ticks[ticks > lim[[i]][1] & ticks < lim[[i]][2]];
ticksMinor <-
ticksMinor[ticksMinor > lim[[i]][1] & ticksMinor < lim[[i]][2]];
if (!tckLab[i]) {
tickLabels <- FALSE;
} else {
tickLabels <- as.character(ticks);
}
# plot the axis
# major ticks
axis(side = i, at = ticks, labels = tickLabels, tck = tck[i],
srt = rotate[[i]]);
# minor ticks
axis(side = i, at = ticksMinor, labels = FALSE, tck = tckMinor[i]);
}
}
invisible(NULL);
}
##~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~.addAxis
## .addAxis2----------------------------2013-03-13
## Modified from the function `.addAxis` to
## add axes to any side of the plot.
## Note: temporary until incorporation and testing.
## ------------------------------------------NB/RH
.addAxis2 <- function (side=1:2, xlim, ylim, tckLab, tck, tckMinor, ...)
{
tckLab <- rep(tckLab, length.out = 4)
tck <- rep(tck, length.out = 4)
tckMinor <- rep(tckMinor, length.out = 4)
mai <- par()$mai
# Nick: The following calls to par() are very strange to me as
# repeated calls to `.addAxis` sequently reduced the size of cex
#par(cex = par()$cex * 0.8)
#par(mai = mai)
#par(mex = par()$cex)
lim <- list(xlim,ylim,xlim,ylim)
rotate <- list(0,90,0,90)
for (i in side) {
if (((i %in% c(1,3)) && (par()$xaxt != "n")) || ((i %in% c(2,4)) && (par()$yaxt != "n"))) {
ticks <- pretty(lim[[i]])
ticksMinor <- pretty(c(0, diff(ticks)[1]))
ticksMinor <- (sort(rep(ticks, length(ticksMinor))) +
rep(ticksMinor, length(ticks)))
ticks <- ticks[ticks > lim[[i]][1] & ticks < lim[[i]][2]]
ticksMinor <- ticksMinor[ticksMinor > lim[[i]][1] &
ticksMinor < lim[[i]][2]]
if (!tckLab[i]) {
tickLabels <- FALSE
}
else {
tickLabels <- as.character(ticks)
}
axis(side = i, at = ticks, labels = tickLabels, tck = tck[i], srt = rotate[[i]], ...)
axis(side = i, at = ticksMinor, labels = FALSE, tck = tckMinor[i], ...)
}
}
invisible(NULL)
}
##~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~.addAxis2
#==============================================================================
.addCorners <- function(polys, ptSummary)
{
xlim <- range(polys$X)
ylim <- range(polys$Y)
corners <- list(tl=c(xlim[1], ylim[2]), tr=c(xlim[2], ylim[2]),
bl=c(xlim[1], ylim[1]), br=c(xlim[2], ylim[1]));
# two "tests" for each corner:
# each list element is the function to apply to X and then Y
tests <- list(tl=c(min, max), tr=c(max, max),
bl=c(min, min), br=c(max, min));
for (c in names(corners)) {
# do we need to add a corner point? check whether a point is already in the
# corner
if (nrow(polys[polys$X == (corners[[c]])[1] &
polys$Y == (corners[[c]])[2], ]) > 0)
next
# find candidate PIDs for the corner;
# candidates touch a boundary on either side of the corner point
polysA <- polys[polys$X == (corners[[c]])[1], ]
if (nrow(polysA) > 0)
# for those polys with an X on the boundary, which ones have a Y
# that is the min/max? get their PIDs
PIDsA <- polysA[polysA$Y == ((tests[[c]])[2])[[1]](polysA$Y), "PID"]
polysB <- polys[polys$Y == (corners[[c]])[2], ]
if (nrow(polysB) > 0)
PIDsB <- polysB[polysB$X == ((tests[[c]])[1])[[1]](polysB$X), "PID"]
if (nrow(polysA) == 0 || nrow(polysB) == 0) {
warning(paste("Unable to close a corner (", c, ").", sep=""));
next
}
# determine candidates
cand <- intersect(PIDsA, PIDsB)
# no candidates (this shouldn't happen...)
if (length(cand) == 0)
stop(paste("Unable to close a corner (", c, ") since no candidates exist.",
sep=""));
# more than one candidate and we need to find the appropriate one
if (length(cand) > 1) {
# compute the distance from each candidate point to the corner
pts <- data.frame (x=ptSummary[cand, "x"], y=ptSummary[cand, "y"])
pts$cornerX <- (corners[[c]])[1];
pts$cornerY <- (corners[[c]])[2];
pts$dist <- sqrt((pts$x - pts$cornerX)^2 + (pts$y - pts$cornerY)^2)
shortest <- which(pts$dist == min(pts$dist))
if (length(shortest) != 1)
stop(paste(
"Unable to determine the appropriate polygon to close corner ", c, ".", sep=""))
cand <- cand[shortest]
}
# add corner point
newPoly <- polys[polys$PID == cand, ]
polys <- polys[polys$PID != cand, ]
xydata <- data.frame(X=c(newPoly$X, corners[[c]][1]),
Y=c(newPoly$Y, corners[[c]][2]));
newPoly <- calcConvexHull(xydata)
newPoly$PID <- cand
polys <- rbind(polys, newPoly);
}
polys <- polys[order(polys$PID), ]
return (polys)
}
#==============================================================================
#.addBubblesLegend----------------------2022-07-05
# Construct legend for the function addBubble.
#--------------------------------------------DC|NB
.addBubblesLegend <- function(radii.leg, usr.xdiff, usr.ydiff,
symbol.zero, symbol.fg, symbol.bg, legend.pos, legend.breaks,
legend.type, legend.title, legend.cex, ...)
{
## ratio of y to x: units-per-inch (Y) / units-per-inch (X)
ratio.y.x = (usr.ydiff / par("pin")[2]) / (usr.xdiff / par("pin")[1])
## calculate the height and width of the legend, which is essential to calculating its position
gap.x <- par("cxy")[1] * legend.cex / 2
gap.y <- par("cxy")[2] * legend.cex / 2
radii.leg.y <- radii.leg * ratio.y.x
leg.tex.w <- strwidth(legend.breaks, units = "user") * legend.cex
title.w = strwidth(legend.title)
max.tex.w <- max(leg.tex.w)
switch(legend.type,
nested = {
## height of the legend is the biggest bubble + title
legend.height <- 2 * max(radii.leg.y) + 3 * gap.y
legend.width <- 2 * max(radii.leg) + gap.x + max.tex.w
},
horiz = {
## height of the legend is the biggest bubble + title + tag
legend.height <- 2 * max(radii.leg.y) + 3 * gap.y
legend.width <- 2 * sum(radii.leg) + (length(legend.breaks) - 1) * gap.x
},
vert = {
legend.height <- 2 * sum(radii.leg.y) + (length(legend.breaks) - 1) * gap.y + 3 * gap.y
legend.width <- 2 * max(radii.leg) + gap.x + max.tex.w
}
)
## reflect an adjustment in X if the title is broader than the legend
if (title.w > legend.width) {
w.adj <- (title.w - legend.width) / 2
} else {
w.adj <- 0
}
## if we already have positions, keep them; otherwise, calculate
## positions given the described corner
#if (class(legend.pos) == "numeric") {
if (inherits(legend.pos, "numeric")) {
legend.loc <- legend.pos
} else {
corners <- c("bottomleft", "bottomright", "topleft", "topright")
if (legend.pos %in% corners) {
legend.loc <- switch(legend.pos,
bottomleft =
c(par("usr")[1] + 0.025 * usr.xdiff + w.adj, par("usr")[3] + 0.025 * usr.ydiff + legend.height),
bottomright =
c(par("usr")[2] - (0.025 * usr.xdiff + legend.width + w.adj), par("usr")[3] + 0.025 * usr.ydiff + legend.height),
topleft =
c(par("usr")[1] + 0.025 * usr.xdiff + w.adj, par("usr")[4] - 0.025 * usr.ydiff),
topright =
c(par("usr")[2] - (0.025 * usr.xdiff + legend.width + w.adj), par("usr")[4] - 0.025 * usr.ydiff)
);
}
}
## the calling function should already have validated legend.type
switch(legend.type,
nested = {
## legend.loc[1] specifies X for *center* of circle; we need
## to shift it right by the largest radius
legend.loc[1] <- legend.loc[1] + max(radii.leg)
## the legend will be drawn from the bottom up; shift
## legend.loc[2] so that it refers to the bottom
legend.loc[2] <- legend.loc[2] - legend.height
r <- rev(radii.leg)
bb <- rev(legend.breaks)
## position of the right edge of the text labels legend
x.text.leg <- legend.loc[1] + max(r) + gap.x + max.tex.w
## draw the circles, lines, and labels
for (i in 1:length(r)) {
symbols(legend.loc[1], legend.loc[2] + r[i] * ratio.y.x,
circles=r[i], inches=FALSE, add=TRUE, bg=symbol.bg[length(r)-i+1], fg=symbol.fg)
lines(c(legend.loc[1], legend.loc[1] + r[1] + gap.x), rep(legend.loc[2] + 2 * r[i] * ratio.y.x, 2))
text(x.text.leg, legend.loc[2] + 2 * r[i] * ratio.y.x, bb[i], adj=c(1, .5), cex=legend.cex)
}
## add the title
x.title.leg <- legend.loc[1] - max(radii.leg) + (legend.width / 2)
text(x.title.leg, legend.loc[2]+legend.height, legend.title, adj=c(0.5,0.5), cex=legend.cex+0.2, col="black")
## set positions for plotting of zero (later)
zlab <- c(x.title.leg, legend.loc[2]+legend.height/4)
},
horiz = {
## legend.loc[2] currently identifies the top of the legend
legend.loc[2] <- legend.loc[2] + max(radii.leg.y) - legend.height
## compute offsets for horizontal spacing
offset <- vector()
for (i in 1:length(radii.leg))
offset[i] <- 2 * sum(radii.leg[1:i]) - radii.leg[i] + (i - 1) * gap.x
## draw circles, labels
symbols(legend.loc[1] + offset, rep(legend.loc[2],length(radii.leg)),
circles = radii.leg, inches = FALSE, bg = symbol.bg, fg = symbol.fg, add = TRUE)
text(legend.loc[1] + offset, legend.loc[2] + radii.leg.y + gap.y,
legend.breaks, adj = c(0.5, 0.5), cex = legend.cex)
## add the title
text(legend.loc[1] + legend.width / 2, legend.loc[2] + legend.height - max(radii.leg.y),
legend.title, adj = c(0.5, 0.5), cex = legend.cex + 0.2, col = "black")
## set positions for plotting of zero (later)
zlab <- c(legend.loc[1], legend.loc[2] - legend.height / 8)
},
vert = {
## part of the calculations have been made above to calculate the height of the legend
if (any(legend.pos == c("bottomleft","topleft")))
legend.loc[1] <- legend.loc[1] + 0.05 * usr.xdiff
## compute offsets for vertical spacing
offset <- vector()
for (i in 1:length(legend.breaks))
offset[i] <- gap.y + 2 * sum(radii.leg.y[1:i]) - radii.leg.y[i] + i * gap.y
## draw circles, labels
symbols(rep(legend.loc[1], length(legend.breaks)), legend.loc[2] - offset,
circles = radii.leg, bg = symbol.bg, fg = symbol.fg, inches = FALSE, add = TRUE)
x.text.leg <- legend.loc[1] + max(radii.leg) + gap.x + max.tex.w
text(rep(x.text.leg, length(legend.breaks)), legend.loc[2] - offset, legend.breaks,
cex = legend.cex, adj = c(1, 0.5), col="black")
## add the title
text(legend.loc[1] + legend.width / 2 - max(radii.leg), legend.loc[2],
legend.title, adj = c(0.5, 0.5), cex = legend.cex + 0.2, col = "black")
## set positions for plotting of zero (later)
zlab <- c(legend.loc[1] + legend.width / 8, legend.loc[2])
}
)
## plot the zero if need be
if (!is.logical(symbol.zero))
legend(zlab[1], zlab[2], legend = "zero", pch = symbol.zero, xjust = 0,
yjust = 1, bty = "n", cex = 0.8, x.intersp = 0.5)
invisible(legend.loc)
}
#==============================================================================
.addFeature <- function(feature, data, polyProps, isEventData,
cex = NULL, col = NULL, font = NULL, pch = NULL, ...)
# 'feature == "points"': add 'data' with points()
# 'feature == "labels"': add 'data' with text()
# 'data': if is PolyData, relaxed unique PID requirement (for addStipples())
# 'isEventData': if TRUE, look for an EID column in 'data'
# '...' contains arguments for either the 'text()' or 'points()' function
#
# Returns: PolyProps
{
data <- .mat2df(data);
if (isEventData) {
type <- "e";
} else {
type <- "p";
}
# given 'feature', appropriate columns in polyProps
if (feature == "points") {
relevantProps <- c("cex", "col", "pch");
} else {
relevantProps <- c("cex", "col", "font");
}
# validate the polyProps argument
polyProps <- .validatePolyProps(polyProps, parCols = relevantProps);
if (is.character(polyProps))
stop(paste("Invalid PolyData 'polyProps'.\n", polyProps, sep=""));
# at this point, 'data' has:
# (!isEventData): PID, (optional SID), X, Y, and label columns
# (isEventData): EID, X, Y, and label columns
# defaults for 'polyProps'
parValues <- list(cex = par("cex"), col = par("col"),
font = par("font"), pch = par("pch"));
# don't replace existing values with defaults
parValues <- parValues[setdiff(names(parValues), names(polyProps))];
# keep only applicable values
parValues <- parValues[intersect(names(parValues), relevantProps)];
# override defaults with values passed as arguments
if (!is.null(cex)) parValues[["cex"]] <- cex;
if (!is.null(col)) parValues[["col"]] <- col;
if (!is.null(font)) parValues[["font"]] <- font;
if (!is.null(pch)) parValues[["pch"]] <- pch;
# make basic 'polyProps' if necessary
if (is.null(polyProps)) {
if (isEventData) {
polyProps <- data.frame(EID = unique(data$EID));
} else {
polyProps <- data.frame(PID = unique(data$PID));
}
}
# merge SIDs into PolyProps if necessary
if (!isEventData
&& is.element("SID", names(data))
&& !is.element("SID", names(polyProps))) {
# grab SIDs from 'polys' for all PIDs in 'polyProps'
p <- data[is.element(data$PID, unique(polyProps$PID)), c("PID", "SID")];
# by filtering things first (above), we speed up the paste() (below)
# add those SIDs to 'polyProps'
polyProps <- merge(polyProps,
p[!duplicated(paste(p$PID, p$SID)), c("PID", "SID")],
by="PID");
}
# flesh out 'polyProps'
if (length(parValues) > 1)
polyProps <- .addProps(type = type, polyProps = polyProps, parValues);
polyPropsReturn <- polyProps;
# reduce data to IDs found in 'polyProps'
if (isEventData) {
data <- data.frame(data[is.element(data$EID, unique(polyProps$EID)), ]);
} else {
if (is.element("SID", names(polyProps))) {
# 'polyProps' may have an SID field and 'polys' may not
if (!is.element("SID", names(data)))
stop("Since 'polyProps' contains an SID column, 'data' must as well.\n");
data <- data[is.element(paste(data$PID, data$SID),
unique(paste(polyProps$PID, polyProps$SID))), ];
} else {
data <- data[is.element(data$PID, unique(polyProps$PID)), ];
}
}
# flatten relevant columns
propColumns <- intersect(names(polyProps), relevantProps);
# paste the columns together
exprStr <- paste("paste(",
paste(paste('polyProps[, "', propColumns, '"]', sep=""),
collapse=", "),
");", sep="");
polyProps$props <- eval(parse(text=exprStr))
# merge that 'props' column into data
if (isEventData) {
data <- merge(data, polyProps[, c("EID", "props")], by="EID");
} else {
# if SID exists in 'polyProps' it exists in 'data' (and vice-versa)
if (is.element("SID", names(polyProps))) {
data <- merge(data,
polyProps[, c("PID", "SID", "props")], by=c("PID", "SID"));
} else {
data <- merge(data, polyProps[, c("PID", "props")], by="PID");
}
}
# split data on the 'props' column
data <- split(data, data$props);
# keep polyProps for getting at individual properties
polyProps <- polyProps[!duplicated(polyProps$props), ];
# plot each set of properties on its own
for (c in names(data)) {
d <- (data[[c]]);
p <- (as.list(polyProps[polyProps$props == c, propColumns]));
if (feature == "labels") {
text(x = d$X, y = d$Y, labels = as.character(d$label),
cex = p$cex, col = p$col, font = p$font, ...);
} else {
points (x = d$X, y = d$Y,
cex = p$cex, col = p$col, pch = p$pch, ...);
}
}
return (polyPropsReturn);
}
## .addLabels---------------------------2008-08-25
## 'projection' should equal "LL" (for longitude/latitude),
## "UTM" (for Universal Transverse Mercator), or anything else (for X/Y).
## 'main', 'sub', 'xlab', 'ylab' may all be part of '...'
## ------------------------------------------NB/RH
.addLabels <- function(projection = NULL, ...)
{
dots <- list(...);
main=dots$main; sub=dots$sub; xlab=dots$xlab; ylab=dots$ylab;
# set label defaults (if necessary)
if (!is.null(projection) && !is.na(projection) && projection == "UTM") {
if (is.null(xlab)) xlab <- "UTM Easting (km)";
if (is.null(ylab)) ylab <- "UTM Northing (km)";
} else if (!is.null(projection) && !is.na(projection) && projection == "LL") {
if (is.null(xlab)) xlab <- "Longitude (\u00B0)"
if (is.null(ylab)) ylab <- "Latitude (\u00B0)"
} else {
if (is.null(xlab)) xlab <- "X";
if (is.null(ylab)) ylab <- "Y";
}
if (is.null(main)) main <- "";
if (is.null(sub)) sub <- "";
# 'xlab'/'ylab' cannot equal NULL at this point
title (main=main, sub=sub, xlab=xlab, ylab=ylab) # cannot add dots - conflicts with formal arguments
invisible(NULL);
}
##~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~.addLabels
## .addProps-----------------------------2004-06-24
## cycles properties in '...' by PID
## 'type = "e"': EventData
## 'type = "p"': PolyData
## 'polyProps': a data frame
## '...': parameters in the form (col = 1, cex = 1, ...),
## each of which can be a list object
## If '...' contains parameters that already exist in 'polyProps', these
## parameters will overwrite the existing columns in 'polyProps' -- avoid
## this behaviour by not including properties in '...' that already exist
## in 'polyProps'.
## Returns: polyProps
##----------------------------------------------NB
.addProps <- function(type, polyProps, ...)
{
# clean up param -- in case it contains one or more lists
param <- list(...);
newParam <- list();
for (i in 1:length(param)) {
if (is.list(param[[i]])) {
newParam <- c(newParam, param[[i]]);
} else {
newParam <- c(newParam, param[i]);
}
}
param <- newParam;
if (type == "e") {
polyProps$IDX <- polyProps$EID;
} else if (type == "p") {
# IDX is PID rather than paste(PID, SID) because we cycle by PID
polyProps$IDX <- polyProps$PID;
} else {
stop (
"Unknown 'type'. Must be either \"e\" or \"p\".");
}
uIDX <- unique(polyProps$IDX);
for (c in names(param)) {
# add it if it isn't NULL
if (!is.null(param[[c]])) {
# remove column if already exists
if (is.element(c, names(polyProps))) {
# use data.frame() to ensure it isn't reduced to a vector
polyProps <- data.frame(polyProps[, !is.element(names(polyProps), c)]);
# if reduced to 1 column, it loses the name; reset it
if (ncol(polyProps) == 1) {
if (type == "p") names(polyProps) <- "PID";
if (type == "e") names(polyProps) <- "EID";
}
}
# build new structure
newColumn <- data.frame(IDX=uIDX);
newColumn[, c] <- rep(param[[c]], length.out = length(uIDX));
# merge new structure
polyProps <- merge(polyProps, newColumn, by = "IDX");
}
}
# remove IDX column; use data.frame() to ensure it isn't reduced to a vector
polyProps <- data.frame(polyProps[, !is.element(names(polyProps), "IDX")]);
# if reduced to 1 column, it loses the name; reset it
if (ncol(polyProps) == 1) {
if (type == "p") names(polyProps) <- "PID";
if (type == "e") names(polyProps) <- "EID";
}
return (polyProps);
}
#==============================================================================
#==============================================================================
.calcDist <- function(polys)
# Assumes 'polys' contains valid PolySet with 'projection' attribute
# containing LL/UTM/1.
#
# Returns: distance vector (distances between each point)
{
# calculate distance for UTM/1:1
if (!is.null(attr(polys, "projection"))
&& !is.na(attr(polys, "projection"))
&& ((attr(polys, "projection") == "UTM")
|| (attr(polys, "projection") == 1))) {
len <- nrow(polys)
D <- c(sqrt((polys$X[1:(len-1)] - polys$X[2:len])^2
+ (polys$Y[1:(len-1)] - polys$Y[2:len])^2),
0);
}
# calculate distance for LL
else if (!is.null(attr(polys, "projection"))
&& !is.na(attr(polys, "projection"))
&& (attr(polys, "projection") == "LL")) {
# Equatorial radius 6,378.14 km
# Polar radius 6,356.78 km
# Mean radius 6,371.3 km
# Sources:
# http://en.wikipedia.org/wiki/Earth
# http://en.wikipedia.org/wiki/Earth_radius
R <- 6371.3;
# degrees to radians
polys[, c("X", "Y")] <- polys[, c("X", "Y")] * pi / 180.0
len <- nrow(polys)
# Source:
# http://www.census.gov/cgi-bin/geo/gisfaq?Q5.1
# Algorithm originally in pseudocode.
s0 <- 1:(len - 1); # 's' for 's'hift
s1 <- 2:len; # (s0 + 1;)
dlon <- polys$X[s1] - polys$X[s0]; # dlon = lon2 - lon1
dlat <- polys$Y[s1] - polys$Y[s0]; # dlat = lat2 - lat1
cosPolysY <- cos(polys$Y);
# a = (sin(dlat/2))^2 + cos(lat1) * cos(lat2) * (sin(dlon/2))^2
a <- (sin(dlat / 2))^2 + cosPolysY[s0] * cosPolysY[s1] *
(sin(dlon / 2))^2;
# c = 2 * arcsin(min(1, sqrt(a)))
a <- sqrt(a);
a[a > 1] <- 1;
cc <- 2 * asin(a);
# d = R * c
D <- c(R * cc, 0);
}
# unknown projection
else {
stop(paste(
"Invalid projection attribute. Supported projections include \"LL\",",
"\"UTM\", and 1.\n"));
}
return (D);
}
#==============================================================================
.calcOrientation <- function(polys)
# Assumes 'polys' contains a valid PolySet.
#
# Returns:
# data frame (invisible) with 'orientation' column (-1 when
# counter-clockwise; 0 when N/A; +1 when clockwise)
# OR: NULL (invisible) (if no rows in output)
{
inRows <- nrow(polys);
# Memory requirement for output are lower than nrow(polys); in fact, they
# are length(unique(paste(polys$PID, polys$SID))).
# The extra time to compute the real memory requirements makes doing so
# not worth it.
outCapacity <- nrow(polys);
# create the data structure that the C functions expect
# Using $<col> notation seems faster than [, "col"] notation.
if (!is.element("SID", names(polys))) {
inID <- c(polys$PID, integer(length = inRows), polys$POS);
} else {
inID <- c(polys$PID, polys$SID, polys$POS);
}
inXY <- c(polys$X, polys$Y);
# call the C function
results <- .C("calcOrientation",
inID = as.integer(inID),
inXY = as.double(inXY),
inVerts = as.integer(inRows),
outID = integer(2 * outCapacity),
outOrientation = double(outCapacity),
outRows = as.integer(outCapacity),
outStatus = integer(1),
PACKAGE = "PBSmapping");
# note: outRows is set to how much space is allocated -- the C function
# should consider this
if (results$outStatus == 1) {
stop(
"Insufficient physical memory for processing.\n");
}
if (results$outStatus == 2) {
stop(paste(
"Insufficient memory allocated for output. Please upgrade to the latest",
"version of the software, and if that does not fix this problem, please",
"file a bug report.\n",
sep = "\n"));
}
# determine the number of rows in the result
outRows <- as.vector(results$outRows);
# extract the data from the C function results
if (outRows > 0) {
d <- data.frame(PID = results$outID[1:outRows],
SID = results$outID[(outCapacity+1):(outCapacity+outRows)],
orientation = results$outOrientation[1:outRows]);
if (!is.element("SID", names(polys)))
d$SID <- NULL;
invisible(d);
} else {
invisible(NULL);
}
}
#==============================================================================
.checkClipLimits <- function(limits)
{
# Makes sure that the X & Y limits are within the bounds of the GSHHS databases
if (limits[1] > limits[2])
stop("xlim[1] is larger than xlim[2]")
if (limits[3] > limits[4])
stop("ylim[1] is larger than ylim[2]")
if (limits[1] > 360 || limits[2] < -20)
stop("xlim are outside of the range of c(-20,360)")
if (limits[3] > 90 || limits[4] < -90)
stop("ylim are outside of the range of c(-90,90)")
}
#==============================================================================
.checkProjection <- function(projectionPlot, projectionPoly)
{
if (is.null(projectionPlot)) {
projMapStr <- "NULL";
} else {
projMapStr <- as.character(projectionPlot);
}
if (is.null(projectionPoly)) {
projPolyStr <- "NULL";
} else {
projPolyStr <- as.character(projectionPoly);
}
msg <- paste(
"The data's 'projection' attribute (", projPolyStr, ") differs from the\n",
"projection of the plot region (", projMapStr, ").\n", sep="");
if (xor(is.null(projectionPlot), is.null(projectionPoly))) {
warning(msg);
} else if ((!is.null(projectionPlot) && !is.null(projectionPoly)) &&
(xor(is.na(projectionPlot), is.na(projectionPoly)))) {
warning(msg);
} else if (!is.null(projectionPlot) && !is.null(projectionPoly) &&
!is.na(projectionPlot) && !is.na(projectionPoly) &&
(projectionPlot != projectionPoly)) {
warning(msg);
}
}
#==============================================================================
.checkRDeps <- function(caller = "unspecified", requires = NULL)
{
if (is.null(version$language) || (version$language != "R")) {
stop (paste ("
The function '", caller, "' requires several dependencies available only in R.\n",
"Please try again from within R.\n", sep=""));
}
err <- NULL;
for (pkg in requires) {
if (!require(pkg, character.only = TRUE)) {
err <- append (err, pkg);
}
}
if (!is.null (err)) {
err <- paste (err, collapse="', '");
stop (paste ("
The function '", caller, "' requires the package(s) '", err, "'.\n",
"Please install the package(s) and try again.\n", sep=""));
}
}
#==============================================================================
.clip <- function(polys, xlim, ylim, isPolygons, keepExtra)
# Does not validate 'polys'; called directly from addPolys since addPolys()
# already validated data.
#
# Maintains extra attributes.
#
# Returns:
# data frame with 'oldPOS' column
# OR: NULL (if no rows in output)
{
# if 'keepExtra', create PolyData from PolySet; after processing, merge these
# data back into 'polys'
if (keepExtra)
pdata <- extractPolyData(polys);
# save the attributes of the data frame
attrNames <- setdiff(names(attributes(polys)),
c("names", "row.names", "class"));
attrValues <- attributes(polys)[attrNames];
inRows <- nrow(polys);
outCapacity <- as.integer(2 * inRows);
# create the data structures that the C function expects
if (!is.element("SID", names(polys))) {
inID <- c(polys$PID, integer(length = inRows), polys$POS);
} else {
inID <- c(polys$PID, polys$SID, polys$POS);
}
inXY <- c(polys$X, polys$Y);
limits <- c(xlim, ylim);
# call the C function
results <- .C("clip",
inID = as.integer(inID),
inXY = as.double(inXY),
inVerts = as.integer(inRows),
polygons = as.integer(isPolygons),
limits = as.double(limits),
outID = integer(4 * outCapacity),
outXY = double(2 * outCapacity),
outRows = as.integer(outCapacity),
outStatus = integer(1),
PACKAGE = "PBSmapping");
# note: outRows is set to how much space is allocated -- the C function
# should take this into consideration
#browser();return() ##***** need to be alert for orphaned holes
if (results$outStatus == 1) {
stop(
"Insufficient physical memory for processing.\n");
}
if (results$outStatus == 2) {
stop(paste(
"Insufficient memory allocated for output. Please upgrade to the latest",
"version of the software, and if that does not fix this problem, please",
"file a bug report.\n",
sep = "\n"));
}
# determine the number of rows in the result
outRows <- as.vector(results$outRows);
# extract the data from the C function results
if (outRows > 0) {
# after creating a data frame several different ways, the following seemed
# to be quickest
d <- data.frame(PID = results$outID[1:outRows],
SID = results$outID[(outCapacity+1):(outCapacity+outRows)],
POS = results$outID[(2*outCapacity+1):(2*outCapacity+outRows)],
oldPOS = results$outID[(3*outCapacity+1):(3*outCapacity+outRows)],
X = results$outXY[1:outRows],
Y = results$outXY[(outCapacity+1):(outCapacity+outRows)]);
# remove "SID" column if not in input
if (!is.element("SID", names(polys)))
d$SID <- NULL;
# change oldPOSs of -1 to NA to meet specs.
d$oldPOS[d$oldPOS == -1] <- NA;
# if keepExtra, merge the PolyData back in...
# 'merge' loses attributes! Merge before restoring attributes
if (keepExtra)
d <- merge(x = d, y = pdata, all.x = TRUE,
by = intersect(c("PID", "SID"), names(d)));
# restore the attributes
attributes(d) <- c(attributes(d), attrValues);
# a final test, to ensure the PolySet contains more than just
# degenerate edges
if (all(d$X == xlim[1]) || all(d$X == xlim[2]) ||
all(d$Y == ylim[1]) || all(d$Y == ylim[2]))
return(NULL);
return(d);
} else {
return(NULL);
}
}
#==============================================================================
.closestPoint <- function(pts, pt)
# "pts": data frame with columns "X" and "Y"
# "pt": data frame with columns "X" and "Y", and ONE row
#
# Returns a vector of length "pts" where T indicates that the point is
# closest to "pt". Returns several Ts when several points equidistant.
{
pts$Xorig <- pt$X
pts$Yorig <- pt$Y
pts$dist <- (pts$Xorig - pts$X) ^ 2 + (pts$Yorig - pts$Y) ^ 2
return (pts$dist == min(pts$dist))
}
#==============================================================================
.createFastIDdig <- function(polysA, polysB = NULL, cols)
# Determines the maximum number of digits in the second column of a data
# frame. If given two data frames ('polysA' and 'polysB', determines that
# maximum between the two data frames.
#
# Assumptions:
# both 'cols' contain integers
#
# Arguments:
# 'polysA': first data frame
# 'polysB': second data frame, which may be missing one or both 'cols';
# if missing one or more 'cols', it's assumed they'll be derived from
# those in 'polysA'
# 'cols': vector of length 2, listing columns to use
#
# Note:
# If it returns 0, then the first and second columns will not fit into
# a double, and they must be 'paste'd toegether.
#
# Returns:
# number of digits in SECOND column (max. of 'polysA' and 'polysB')
# OR NULL if the columns don't exist in the data frame(s)
{
if ((length(cols) == 2)
&& all(is.element(cols, names(polysA)))) {
digitsL <- floor(log10(max(polysA[[cols[1]]])) + 1);
digitsR <- floor(log10(max(polysA[[cols[2]]])) + 1);
# check 'polysB' as well
if (!is.null(polysB)) {
if (is.element(cols[1], names(polysB)))
digitsL <- max(digitsL, floor(log10(max(polysB[[cols[1]]])) + 1));
if (is.element(cols[2], names(polysB)))
digitsR <- max(digitsR, floor(log10(max(polysB[[cols[2]]])) + 1));
}
# 'double' has 15 digits of precision (decimal), according to my
# 'C Pocket Reference' by Prinz, P. and U. Kirch-Prinz
if ((digitsL + digitsR) <= 15) {
return (digitsR);
} else {
return (0);
}
} else {
return (NULL);
}
}
#==============================================================================
.createGridIDs <- function(d, addSID, byrow)
# Create IDs for a grid according to the addSID and byrow arguments.
#
# Arguments:
# 'd': PolySet (grid) created using addSID = T and byrow = T
# 'addSID': if TRUE, include an SID column
# 'byrow': if TRUE, increment PID along X
{
if (addSID && !byrow) {
# swap
tmp <- d$PID;
d$PID <- d$SID;
d$SID <- tmp;
} else if (!addSID && byrow) {
d$PID <- (d$SID - 1) * (length(unique(d$X)) - 1) + d$PID;
d$SID <- NULL;
} else if (!addSID && !byrow) {
d$PID <- (d$PID - 1) * (length(unique(d$Y)) - 1) + d$SID;
d$SID <- NULL;
}
return (d);
}
#==============================================================================
.createIDs <- function(x, cols, fastIDdig = NULL)
# Creates an IDs (or IDX) column from its input.
#
# Arguments:
# 'x': data frame with one or more columns
# 'cols': columns to use when creating the index; OK if individual
# columns are missing
# 'fastIDdig': (optional) maximum number of digits in the second column, if
# only two integer columns, and the maximum number of digits between them
# is less than 15; often the output from '.createFastIDdig'; a user would
# want to pass in 'fastIDdig' when creating (some) matching indices for two
# different data frames
#
# Note:
# If 'fastIDdig' equals NULL, executes '.createFastIDdig' to create it (but
# only if two columns).
#
# Returns:
# index column if everything OK
# NULL if error
{
# use 'is.element' instead of 'intersect' because 'intersect' gives no
# guarantee of order (as far as I can see)
presentCols <- cols[is.element(cols, names(x))]
if (length(presentCols) == 1) {
return (x[[presentCols]]);
} else if (length(presentCols) == 2) {
# if a fastIDdig wasn't passed in, try to create one
if (is.null(fastIDdig)) {
fastIDdig <- .createFastIDdig(polysA=x, polysB=NULL, cols=presentCols);
}
# if called the function above, 'fastIDdig' won't equal NULL (but might
# equal 0)
# if the user passed it in, it should only equal NULL if there is
# only one ID column, in which case we shouldn't be in this branch
# of the 'if'
if (fastIDdig > 0) {
return (as.double(x[[presentCols[1]]]
+ (x[[presentCols[2]]] / 10^fastIDdig)));
} else {
return (paste(x[[presentCols[1]]], x[[presentCols[2]]], sep = "-"));
}
} else if (length(presentCols) > 2) {
exprStr <- paste("paste(",
paste(paste("x$", presentCols, sep=""),
collapse=", "),
");", sep="");
return (eval(parse(text=exprStr)));
}
return (NULL);
}
#==============================================================================
.expandEdges <- function(polys, pts, xlim, ylim)
# expects all arguments to be specified
{
polyRange <- c(range(polys$X), range(polys$Y))
ptsRange <- c(range(pts$X), range(pts$Y))
# set "toFix" based on points outside the PolySet
toFixPts <- c(ptsRange[1] < polyRange[1],
ptsRange[2] > polyRange[2],
ptsRange[3] < polyRange[3],
ptsRange[4] > polyRange[4])
# update "toFix" to extend PolySet to limits outside of its range
toFixLim <- c(signif(xlim[1], 5) < signif(polyRange[1], 5),
signif(xlim[2], 5) > signif(polyRange[2], 5),
signif(ylim[1], 5) < signif(polyRange[3], 5),
signif(ylim[2], 5) > signif(polyRange[4], 5))
# the use of "signif" reduces the number of unnecessary edges in
# "toFixLim" due to rouding error; digits = 6 seemed insufficient
# to compensate for rounding error
toFix <- toFixPts | toFixLim
# if nothing to fix, return
if (!any(toFix))
return (polys)
# side: 1 = left, 2 = right, 3 = bottom, 4 = top
for (side in which(toFix)) {
# a point that causes expansion
if (side == 1) {
PID <- which(pts$X < polyRange[side])
} else if (side == 2) {
PID <- which(pts$X > polyRange[side])
} else if (side == 3) {
PID <- which(pts$Y < polyRange[side])
} else if (side == 4) {
PID <- which(pts$Y > polyRange[side])
} else {
stop ("Internal error: unrecognized value of \"size\" in point check.")
}
# no points caused the expansion, so it must be a limit that causes it
if (length(PID) == 0) {
# must determine PID using another method
if (side == 1 || side == 2) {
PID <- which(.closestPoint(pts, data.frame(X=xlim[side], Y=mean(ylim))))
} else if (side == 3 || side == 4) {
PID <- which(.closestPoint(pts, data.frame(X=mean(xlim), Y=ylim[side - 2])))
} else {
stop ("Internal error: unrecognized value of \"size\" in limit check.")
}
}
# sanity check
if (length(PID) != 1)
stop ("Internal error: unable to determine appropriate PID for expansion.")
newPoly <- polys[polys$PID == PID, ]
polys <- polys[polys$PID != PID, ]
if (side == 1 || side == 2) {
newXY <- data.frame(X=c(newPoly$X, rep(xlim[side], 2)),
Y=c(newPoly$Y, range(newPoly$Y)))
} else if (side == 3 || side == 4) {
newXY <- data.frame(X=c(newPoly$X, range(newPoly$X)),
Y=c(newPoly$Y, rep(ylim[side - 2], 2)))
} else {
stop ("Internal error: unrecognized value of \"size\" in data frame setup.")
}
newPoly <- calcConvexHull(newXY)
newPoly$PID <- PID
polys <- rbind(polys, newPoly)
}
polys <- polys[order(polys$PID), ]
return (polys);
}
#==============================================================================
.fixGSHHSWorld <- function (world) {
# store desired limits
xlim <- range(world$X)
ylim <- range(world$Y)
ylim[1] <- -90
# determine PID of Antarctica: we'll use it to extract the current Antarctica,
# which we'll grow west/east, then clip, then merge with the other polygons
event <- data.frame(EID = 1, X = 85, Y = -72)
event <- findPolys(event, world)
pid <- event$PID[1]
# extract current Antarctica
curAnt <- world[is.element(world$PID,pid), ]
# ensure the points are ordered from left to right
if (curAnt$X[1] > curAnt$X[nrow(curAnt)])
curAnt <- curAnt[nrow(curAnt):1, ]
# create a copy to the left and a copy to the right; add corners
left <- curAnt[c(1,1:nrow(curAnt)), ]
left$X <- left$X - 360
left$Y[1] <- -90
right <- curAnt[c(1:nrow(curAnt), nrow(curAnt)), ]
right$X <- right$X + 360
right$Y[nrow(right)] <- -90
# merge to create a new, very wide Antarctica
curAnt <- rbind(left, curAnt, right)
curAnt$POS <- 1:nrow(curAnt)
# clip
curAnt <- clipPolys(curAnt, xlim, ylim)
curAnt$oldPOS <- NULL
# merge into the existing world
world <- rbind(world[world$PID < pid, ], curAnt, world[world$PID > pid, ])
row.names(world) <- 1:length(row.names(world))
invisible(world)
}
#=============================================================================
.getBasename <- function (fn, ext)
# If appropriate, remove the extension from 'fn' to obtain the shapefile
# name without the extension. When testing for file existance, use the
# extension 'ext'.
{
# if appending .shp does not give a valid file...
if (!file.exists (paste(fn, ".", ext, sep=""))) {
# attempt to remove extension and try again
fn <- sub ("\\..{3}$", "", fn);
if (!file.exists (paste(fn, ".", ext, sep=""))) {
stop (paste ("Cannot find the file \"", fn, ".", ext, "\".", sep=""));
}
}
return (fn);
}
#=============================================================================
.getGridPars <- function (polys, fullValidation = TRUE)
{
res <- list();
res$x <- sort(unique(polys$X));
res$y <- sort(unique(polys$Y));
lenx <- length(res$x);
leny <- length(res$y);
if ((lenx < 2) || (leny < 2)) {
return (NULL);
}
# determine "addSID"
if (is.element("SID", names(polys))) {
res$addSID <- TRUE;
} else {
res$addSID <- FALSE;
}
# determine "byrow"
if (lenx == 2 && leny == 2) {
# special case: only one polygon; byrow does not matter
res$byrow <- TRUE;
} else if (lenx == 2) {
# special case: the second polygon in polys is _above_ the
# first
if (polys[1, "PID"] + 1 == polys[5, "PID"]) {
res$byrow <- FALSE;
} else {
res$byrow <- TRUE;
}
} else {
# the second polygon in polys is to the _right_ of the
# first
if (polys[1, "PID"] + 1 == polys[5, "PID"]) {
res$byrow <- TRUE;
} else {
res$byrow <- FALSE;
}
}
res$projection <- attr(polys, "projection");
res$zone <- attr(polys, "zone");
if (fullValidation) {
t <- makeGrid(x = res$x, y = res$y, byrow = res$byrow,
addSID = res$addSID, projection = res$projection,
zone = res$zone);
if (is.character (all.equal (polys, t))) {
return (NULL);
}
}
return (res);
}
#==============================================================================
.initPlotRegion <- function(projection, xlim, ylim, plt)
# Initialize the plot region, accounting for the aspect ratio.
#
# 'projection': projection
# "UTM" => UTM <-- as per specs.
# "LL" => Geographic <-- as per specs.
# NA => none <-- for plotLines, plotPolys
# numeric => specifies aspect ratio
# 'xlim': x-limits of the plot
# 'ylim': y-limits of the plot
# 'plt': requested plot size
#
# Returns: NULL (invisible)
#
# Notes:
# In a new plot, mar/mai are 'in-sync' with pin/plt, meaning that changing
# mar/mai causes pin/plt to change. Manually setting pin/plt causes mar/mai
# to fall 'out-of-sync', meaning that changing mar/mai will have no effect
# on pin/plt. In this function, do not cause pin/plt to fall out-of-sync
# (i.e., never directly set pin/plt). However, if they are already
# out-of-sync, honour the pin setting (which is consistent with what R/S-PLUS
# does).
#
# Summary:
# Read plt; whether in-sync or out-of-sync, this value reflects the
# maximum plot region. Calculate desired plot region, assume in-sync, and
# set new region indirectly via mar/mai. If new region doesn't match the
# expected region, it is out-of-sync, so set new region directly via pin/plt.
#
# If the user allows pin/plt to fall out-of-sync and sets plt to NULL, the
# function will still produce the correct aspect ratio, but may not honour
# the mai/mar par() parameters.
{
# validate 'xlim' and 'ylim'
if ((missing(xlim) || is.null(xlim)) || !is.vector(xlim) ||
(length(xlim) != 2) ||
(missing(ylim) || is.null(ylim)) || !is.vector(ylim) ||
(length(ylim) != 2)) {
stop(
"xlim and/or ylim is missing or invalid.\n");
}
if (diff(ylim) == 0 || diff(xlim) == 0) {
stop(
"xlim/ylim must specify a region with area greater than 0.\n");
}
# since all high-level plot functions call this function, we can
# advance to the next frame here, unless par(new=T) is specified.
if(!par()$new)
frame();
# if 'plt' is NULL, the user doesn't want to use the default, so read its
# current value from 'par' parameters
if (is.null(plt))
plt <- par("plt");
# 'projection' cannot equal NULL
if (is.null(projection))
stop(
"'projection' argument must not equal NULL.\n");
# if projection isn't NA, adjust for the aspect ratio; in this section
# we must indirectly set plt/pin via mai/mar if possible; if that doesn't
# work, we must set them directly
if (!is.na(projection))
{
# use an additional scaling factor if LL
xyRatio <- ifelse(projection == "LL", cos((mean(ylim) * pi) / 180), 1);
# determine plot region aspect ratio; calculate indirectly using
# fin and plt
aspPlotRegion <-
(par()$fin[1]*diff(plt[1:2])) / (par()$fin[2]*diff(plt[3:4]));
if (is.infinite(aspPlotRegion))
stop("Plot region must have an area greater than 0.\n");
# determine desired aspect ratio; xyRatio != 1 only when proj == "LL"
aspPolySet <- diff(xlim) / diff(ylim) * xyRatio;
if (is.numeric(projection))
aspPolySet <- aspPolySet / projection;
if (is.infinite(aspPolySet) || (aspPolySet == 0))
stop(paste(
"Either 'projection' is 0 or 'xlim'/'ylim' specify a region with an area",
"of 0.\n", sep="\n"));
# adjust 'plt' to honour the aspect ratio
if (aspPlotRegion < aspPolySet) {
# shrink the Y direction
pinX <- par()$fin[1]*diff(plt[1:2]);
pinY <- pinX * aspPolySet^-1;
toMove <- (diff(plt[3:4]) - pinY/par()$fin[2]) / 2;
plt[3:4] <- plt[3:4] + c(toMove, -toMove);
} else if (aspPlotRegion > aspPolySet) {
# shrink the X direction
pinY <- par()$fin[2]*diff(plt[3:4]);
pinX <- pinY * aspPolySet;
toMove <- (diff(plt[1:2]) - pinX/par()$fin[1]) / 2;
plt[1:2] <- plt[1:2] + c(toMove, -toMove);
}
# set the plot region indirectly via mai/mar
# mai: (bottom, left, top, right)
# fin: (width, height)
# plt: (left, right, bottom, top)
# If out-of-sync, may issue a warning re: plot specified in inches
# too large.
par(mai=c(par()$fin[2]*plt[3], par()$fin[1]*plt[1],
par()$fin[2]*(1-plt[4]), par()$fin[1]*(1-plt[2])));
# verify that it set OK; if not, set it directly via plt
parPlt <- signif(as.double(par()$plt), digits=5);
locPlt <- signif(as.double(plt), digits=5);
if (any(parPlt != locPlt))
par(plt = plt);
}
# if projection is NA, continue at this point with no changes to
# the plot region
# must set 'usr' again because changing the other par() parameters changes it
par(usr = c(xlim, ylim));
invisible(NULL);
}
#==============================================================================
.insertNAs <- function(polys, idx) {
# extract the polys of interest
sel.polys <- polys[is.element(names(polys), idx)]
# NOTE: the above code is _much_ faster in S-PLUS than an alternative...
# sel.polys <- polys[as.character(idx)]
# ... which I tried first
if (is.null(unlist(sel.polys)))
return (NA)
# lapply's are slower in S-PLUS than R; the next line is a bottle-neck
# in S-PLUS (while having minimal impact on R runtime)
lenPolys <- lapply(sel.polys, "length")
nPolys <- length(lenPolys)
# create a vector T, F, ..., T
TFT <- rep(c(TRUE, FALSE), length.out=(nPolys * 2) - 1)
# create a vector with the number of vertices/polygon, with each separated
# by a 1
reps <- rep(1, len=(nPolys * 2) - 1)
reps[TFT] <- unlist(lenPolys)
# create a vector that can select all of the NAs that seperate polygons
# in the "new.polys" vector
NAs <- rep(!TFT, times=reps)
new.polys <- vector(length=(length(unlist(sel.polys)) + nPolys - 1))
new.polys[NAs] <- NA
new.polys[!NAs] <- unlist(sel.polys)
return (new.polys)
}
## .is.in-------------------------------2018-09-06
## Check to see if a set of vertices occurs in a
## polygon (based on sp::point.in.polygon C code)
## ------------------------------------------NB|RH
.is.in <- function(events, polys, use.names=TRUE)
{
if (!is.PolySet(polys))
stop("Supply a PolySet to argument 'polys'")
if (!all(c(is.element(c("X","Y"),colnames(events)),is.element(c("X","Y"),colnames(polys)))))
stop("Objects 'events' and 'polys' must have column names 'X' and 'Y'")
## Code lifted from NB's code in `findPolys'
## Create the data structure that the C function expects
if (!is.EventData(events))
events = as.EventData(data.frame(EID=1:nrow(events), events), projection=attributes(polys)$projection)
inEvents = nrow(events)
inEventsID = events$EID
inEventsXY = c(events$X, events$Y)
inPolys = nrow(polys)
if (!is.element("SID", names(polys))) {
inPolysID = c(polys$PID, integer(length=inPolys), polys$POS)
} else {
inPolysID = c(polys$PID, polys$SID, polys$POS)
}
inPolysXY = c(polys$X, polys$Y)
## The maximum number of rows in the results occurs when each event
## lies in all of the polygons; we can easily calculate this case
## but it often results in too much allocated memory
## Instead, simply use an argument that the user can tune (this makes
## the behavior more consistent with joinPolys, too)
outCapacity <- nrow(events)
## Call the C function
results <- .C("findPolys",
inEventsID = as.integer(inEventsID),
inEventsXY = as.double(inEventsXY),
inEvents = as.integer(inEvents),
inPolysID = as.integer(inPolysID),
inPolysXY = as.double(inPolysXY),
inPolys = as.integer(inPolys),
outID = integer(4 * outCapacity),
outRows = as.integer(outCapacity),
outStatus = integer(1),
PACKAGE = "PBSmapping")
## Note: outRows is set to how much space is allocated -- the C function should consider this
## Determine the number of rows in the result
outRows <- as.vector(results$outRows)
if (outRows == 0) {
e.out = events
e.in = events[0,]
e.bdry = 0
} else {
## Extract the data from the C function results
d <- data.frame(EID = results$outID[1:outRows],
PID = results$outID[(outCapacity+1):(outCapacity+outRows)],
SID = results$outID[(2*outCapacity+1):(2*outCapacity+outRows)],
Bdry = results$outID[(3*outCapacity+1):(3*outCapacity+outRows)])
e.in = events[is.element(events$EID,d$EID),]
e.out = events[!is.element(events$EID,d$EID),]
e.bdry = d$Bdry
}
out= list()
out[["e.in"]] = e.in
out[["e.out"]] = e.out
out[["all.in"]] = nrow(e.in)==nrow(events)
out[["all.out"]] = nrow(e.out)==nrow(events)
out[["all.bdry"]] = all(e.bdry==1)
return(out)
#browser();return()
}
##~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~.is.in
### .is.in.defunct-------------------------------2018-08-24
### Check to see if a set of vertices occurs in a
### polygon (based on sp::point.in.polygon C code)
#.is.in.defunct =function(events, polys, use.names=TRUE) {
# if (!all(c(is.element(c("X","Y"),colnames(events)),is.element(c("X","Y"),colnames(polys)))))
# stop("Objects 'events' and 'polys' must have column names 'X' and 'Y'")
# .checkRDeps(".is.in", c("sp"))
# eval(parse(text="in.out = .Call(\"R_point_in_polygon_sp\", as.numeric(events[,\"X\"]), as.numeric(events[,\"Y\"]), as.numeric(polys[,\"X\"]), as.numeric(polys[,\"Y\"]), PACKAGE = \"sp\")"))
# if (use.names) {
# if(is.PolySet(events) || is.PolyData(events))
# names(in.out) = .createIDs(events, intersect(c("PID","SID"), colnames(events)))
# if(is.EventData(events))
# names(in.out) = events[,"EID"]
# }
# out = list()
# out[["in.out"]] = in.out
# out[["all.in"]] = all(as.logical(in.out))
# return(out)
##browser();return()
#}
###~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~.is.in.defunct
#==============================================================================
.mat2df <- function(data) {
# Convert matrices to data.frames, preserving certain attributes.
# Attributes to ignore: dim, dimnames, class
pbsClass <- intersect(attributes(data)$class, c("EventData", "LocationSet",
"PolyData", "PolySet"));
attrNames <- names(attributes(data));
if(all(is.element(attrNames, c("dim", "dimnames", "class")) == TRUE)) {
addValues <- NULL;
} else {
addNames <- setdiff(attrNames, c("dim", "dimnames", "class"));
addValues <- attributes(data)[addNames];
}
data <- data.frame(unclass(data));
if(!is.null(addValues)) {
attributes(data) <- c(attributes(data), addValues);
}
if(length(pbsClass) > 0) {
attr(data, "class") <- c(pbsClass, "data.frame");
}
return(data);
}
#==============================================================================
.rollupPolys <- function(polys, rollupMode, exteriorCCW, closedPolys,
addRetrace)
# rollupMode: method for rolling up the PolySet; essentially controls when
# to introduce a new PID (or PID/SID)
# 1 = roll-up to the PID level (only PIDs in the result)
# 2 = roll-up to the outer contour level (only outer contours in the
# result)
# 3 = do not roll-up
#
# exteriorCCW: modify vertices orientation (CW/CCW)?
# -1 = don't modify
# 0 = exterior should be CW
# +1 = exterior should be CCW
#
# closedPolys: whether the last and first vertices should be the same
# -1 = don't modify
# 0 = ensure polygons do not close
# +1 = close the polygons
#
# addRetrace: determines whether it adds retrace lines to the first vertex
# of the parent after outputting a child
# 0 = don't add
# 1 = add
#
# Note: does not validate the PolySet
#
# Returns: PolySet (invisible) or NULL (invisible)
{
# save the attributes for the data frame (.validatePolySet returns a data
# frame)
attrNames <- setdiff(names(attributes(polys)),
c("names", "row.names", "class"));
attrValues <- attributes(polys)[attrNames];
inRows <- nrow(polys);
outCapacity <- 2 * inRows;
# create the data structure that the C function expects
if (!is.element("SID", names(polys))) {
inID <- c(polys$PID, integer(length = inRows));
} else {
inID <- c(polys$PID, polys$SID);
}
inXY <- c(polys$X, polys$Y);
results <- .C("rollupPolys",
inID = as.integer(inID),
inPOS = as.double(polys$POS),
inXY = as.double(inXY),
inVerts = as.integer(inRows),
outID = integer(3 * outCapacity),
outXY = double(2 * outCapacity),
outRows = as.integer(outCapacity),
rollupMode = as.integer(rollupMode),
exteriorCCW = as.integer(exteriorCCW),
closedPolys = as.integer(closedPolys),
addRetrace = as.integer(addRetrace),
outStatus = integer(1),
PACKAGE = "PBSmapping");
# note: outRows is set to how much space is allocated -- the C function
# should take this into consideration
if (results$outStatus == 1) {
stop(
"Insufficient physical memory for processing.\n");
}
if (results$outStatus == 2) {
stop(paste(
"Insufficient memory allocated for output. Please upgrade to the latest",
"version of the software, and if that does not fix this problem, please",
"file a bug report.\n",
sep = "\n"));
}
if (results$outStatus == 3) {
stop(paste(
"Unable to rollup the polygons, as one or more children did not have a",
"parent.\n"));
}
# determine the number of rows in the result
outRows <- as.vector(results$outRows);
# extract the data from the C function results
if (outRows > 0) {
d <- data.frame(PID = results$outID[1:outRows],
SID = results$outID[(outCapacity+1):(outCapacity+outRows)],
POS = results$outID[(2*outCapacity+1):(2*outCapacity+outRows)],
X = results$outXY[1:outRows],
Y = results$outXY[(outCapacity+1):(outCapacity+outRows)]);
if (!is.element("SID", names(polys)) || rollupMode == 1)
d$SID <- NULL;
# restore the attributes
attributes(d) <- c(attributes(d), attrValues);
invisible(d);
} else {
invisible(NULL);
}
}
#=====================================================================2018-06-05
# Default for 'projection' must be logical type.
# For 'plotPoints()', 'cex'and 'pch' are part of '...' -- and they may equal
# NULL; when they equal NULL, don't try adding them to par()!
.plotMaps <- function(polys, xlim, ylim, projection, plt, polyProps,
border, lty, col, colHoles, density, angle, bg,
axes, tckLab, tck, tckMinor, isType, ...)
{
# - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
# PART 1: deal with the '...' business; backup original 'par' values, and set
# values specified in '...' in 'par'
# DO NOT backup 'new' (and some others)
legalNames <- c("adj", "ann", "ask", "bg", "bty", "cex", "cex.axis",
"cex.lab", "cex.main", "cex.sub", "col", "col.axis",
"col.lab", "col.main", "col.sub", "crt", "csi", "err",
"exp", "fg", "font", "font.axis", "font.lab",
"font.main", "font.sub", "lab", "las", "lty",
"lwd", "mgp", "mkh", "pch", "smo", "srt", "tck",
"tcl", "tmag", "type", "xaxp", "xaxs", "xaxt", "xpd",
"yaxp", "yaxs", "yaxt");
if (!is.null(version$language) && (version$language == "R")) {
legalNames <- setdiff(legalNames, "csi"); # read-only in R
}
legalNames <- intersect(legalNames, names(par()));
backupPar <- par(legalNames);
on.exit(par(backupPar));
dots <- list(...);
extraArgs <- setdiff(names(dots), legalNames);
extraArgs <- setdiff(extraArgs, c("main", "sub", "type", "xlab", "ylab"));
if (length(extraArgs) > 0) {
warning(paste(
"Ignored unrecognized argument '", paste(extraArgs, collapse="', '"), "'.\n",
sep = ""));
}
# special case: since they aren't explicit parameters, they end up being
# part of '...'; but they may equal NULL (invalid values); remove them, since
# they are explicitly passed into addPoints(), where they become part
# of 'polyProps'
if (isType == "points")
legalNames <- setdiff(legalNames, c("cex", "pch"));
par(dots[intersect(names(dots), legalNames)]);
# - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
# PART 2: set 'projection' and 'labelprojection' to the appropriate values
# Since we validate here, it's unnecessary to validate when this function
# calls 'addLines'/'addPolys'
# although 'projection' cannot equal NA from caller, we may set it
# to NA within this function (to mean unprojected)
if (is.null(projection) || is.na(projection))
projection <- FALSE;
if (is.logical(projection)) {
if (projection) {
if (!is.null(attr(polys, "projection"))) {
projection <- attr(polys, "projection");
} else {
projection <- 1;
warning(
"'projection' set to 1:1 since unspecified 'projection' argument/attribute.\n");
}
} else {
projection <- NA;
}
} else {
# user specified 'projection': should be numeric, "LL", "UTM", and should
# not conflict with 'polys'
if (is.numeric(projection) || is.element(projection, c("LL", "UTM"))) {
if (!is.null(attr(polys, "projection")) &&
!is.element(projection, attr(polys, "projection"))) {
projection <- attr(polys, "projection");
warning(
"'projection' argument overwritten with PolySet's 'projection' attribute.\n");
}
} else {
stop(paste(
"Either omit 'projection' argument or set it to a numeric value, \"LL\", or",
"\"UTM\".\n",
sep="\n"));
}
}
# prevent duplicate warnings on call to 'addLines'/'addPolys'
if (!is.null(polys))
attr(polys, "projection") <- projection;
labelProjection <- projection;
# 'projection' type no longer 'logical'; now
# NA: no projection
# numeric: specified projection
# "LL": longitude-latitude
# "UTM": UTM
# - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
# PART 3: get on with the routine
# special case of NULL PolySet
if (is.null(polys)) {
if (is.null(xlim) || is.null(ylim) || is.null(projection)) {
stop(
"To plot a NULL PolySet, pass 'xlim', 'ylim', and 'projection' arguments.\n");
}
}
else {
# as far as this function is concerned, 'polys' must have an X/Y, and
# that's it; do a basic check
polys <- .validateXYData(polys);
if (is.character(polys)) stop(paste("Invalid PolySet.\n", polys, sep=""));
}
# detect limits if necessary
if (is.null(xlim)) xlim <- range(polys$X);
if (is.null(ylim)) ylim <- range(polys$Y);
# check 'dots' for 'type'
if (is.element("type", names(dots))) {
if (dots$type == "n") {
# set 'polys' to NULL since we've decided that type == "n" is equivalent
# to passing a NULL PolySet
polys <- NULL;
} else {
stop(
"Either omit 'type' argument or set it to \"n\".\n");
}
}
# set 'col' parameter appropriately for adding labels/etc.
if (!is.null(polys)) {
par(col = 1);
} else {
if (length(col) > 1) {
stop(paste(
"Either omit 'col' argument or set it to a single-element vector when 'polys'",
"equals NULL or 'type = \"n\"'.\n",
sep = "\n"));
} else if (!is.null(col)) {
par(col = col);
}
}
# save settings in 'options'
options(map.xlim = xlim);
options(map.ylim = ylim);
options(map.projection = projection);
# create plot region
.initPlotRegion(projection=projection, xlim=xlim, ylim=ylim, plt=plt);
## plot background colour
## note (180605): RH had to specify xpd=TRUE to get `polygon' to honour the pin regions on Windows devga device.
## xpd: a logical value or NA.
## if FALSE, all plotting is clipped to the plot region (appears to be buggy),
## if TRUE, all plotting is clipped to the figure region (adopted herein), and
## if NA, all plotting is clipped to the device region (can also use this but overkill).
if (!is.null(bg))
polygon(x=xlim[c(1,2,2,1)], y=ylim[c(1,1,2,2)], col=bg, border=0, xpd=TRUE);
# plot PolySet 'polys'
if (!is.null(polys)) {
if (isType == "polygons") {
# add polygons
ret <- addPolys(polys, xlim = xlim, ylim = ylim, polyProps = polyProps,
border = border, lty = lty, col = col, colHoles = colHoles,
density = density, angle = angle);
} else if (isType == "lines") {
# add lines
ret <- addLines(polys, xlim = xlim, ylim = ylim, polyProps = polyProps,
lty = lty, col = col);
} else if (isType == "points") {
cex <- list(...)$cex; # these were passed in '...'
pch <- list(...)$pch;
ret <- addPoints(polys, xlim = xlim, ylim = ylim, polyProps = polyProps,
cex = cex, col = col, pch = pch);
} else {
stop(
"Unrecognized 'isType'.\n");
}
} else {
ret <- NULL;
}
if (axes) {
.addAxis(xlim = xlim, ylim = ylim, tckLab = tckLab, tck = tck,
tckMinor = tckMinor, ...);
}
else {
options(map.xline = 1);
options(map.yline = 1);
}
# labels must go after axis
.addLabels(projection = labelProjection, ...);
if (axes) {
# since R won't plot outside the figure region, we don't want the
# box to have a thicker line width
box();
}
invisible(ret);
}
#==============================================================================
# .preparePolyProps:
# Performs at least the following tasks:
# 1) creates 'polyProps' if it equals NULL
# 2) adds an SID column to 'polyProps' if one exists in 'polys'
# 3) removes from 'polyProps' any PIDs that do not exist in 'polys'
.preparePolyProps <- function(polysPID, polysSID, polyProps)
{
# make 'polyProps' if necessary
if (is.null(polyProps)) {
polyProps <- data.frame(PID = unique(polysPID));
} else {
polyProps <- .validatePolyData(polyProps);
if (is.character(polyProps))
stop(paste("Invalid PolyData 'polyProps'.\n", polyProps, sep=""));
}
# if SIDs in 'polys' but not 'polyProps', add SIDs from 'polys' into
# 'polyProps'
if (!is.null(polysSID) && !is.element("SID", names(polyProps))) {
# identify unique PIDs/SIDs within "polys"
p <- data.frame(PID=polysPID, SID=polysSID);
p <- p[!duplicated(.createIDs(p, cols = c("PID", "SID"))), ]
# only keep PIDs that appear in 'polyProps'
p <- p[is.element(p$PID, unique(polyProps$PID)), ];
# add the SIDs to 'polyProps'
polyProps <- merge(polyProps, p, by="PID");
}
return (polyProps);
}
## .validateData------------------------2018-09-10
## An element of noNACols and keyCols will only be
## used if it exists in the data. To ensure it
## exists in the data, make it a requiredCol.
## ------------------------------------------NB/RH
.validateData <- function(data, className,
requiredCols=NULL, requiredAttr=NULL, noFactorCols=NULL,
noNACols=NULL, keyCols=NULL, numericCols=NULL)
{
fnam = as.character(substitute(data))
## convert matrix to data frame
if (is.matrix(data)) {
data <- .mat2df(data);
}
if (is.data.frame(data) && (nrow(data) > 0)) {
## validate optional class name
if (!is.null(className) && (class(data)[1] != "data.frame")) {
if (class(data)[1] != className) {
return(paste("Unexpected class (", class(data)[1], ").\n", sep=""));
}
}
## ensure all the required columns exist in the PolySet
if (!is.null(requiredCols) && !all(is.element(requiredCols, names(data)))) {
return(paste("One or more of the required columns is missing.\n",
"Required columns: ", paste(requiredCols, collapse = ", "), ".\n", sep=""));
}
## ensure all the required attributes exists in the PolySet
if (!is.null(requiredAttr) && !all(is.element(requiredAttr, names(attributes(data))))) {
return(paste("One or more of the required attributes is missing.\n",
"Required attributes: ", paste(requiredAttr, collapse = ", "), ".\n", sep=""));
}
## check for NAs
presentCols <- intersect(noNACols, names(data));
if (length(presentCols) > 0) {
# build an expression
exprStr <- paste(paste("any(is.na(data$", presentCols, "))", sep=""), collapse=" || ");
if (eval(parse(text=exprStr))) {
return(paste("One or more columns (where NAs are not allowed) contains NAs.\n",
"Columns that cannot contain NAs: ", paste(presentCols, collapse = ", "),
".\n", sep=""));
}
}
## check for factors
presentCols <- intersect(noFactorCols, names(data));
if (length(presentCols) > 0) {
# build an expression
exprStr <- paste(paste("is.factor(data$", presentCols, ")", sep=""), collapse=" || ");
if (eval(parse(text=exprStr))) {
return(paste("One or more columns contains factors where they are not allowed.\n",
"Columns that cannot contain factors: ", paste(presentCols, collapse = ", "),
".\n", sep=""));
}
}
## check for uniqueness of the keys
presentCols <- intersect(keyCols, names(data));
if (length(presentCols) > 0) {
if (length(presentCols) == 1) {
keys <- data[[presentCols]];
} else if ((length(presentCols) == 2) && ((length(intersect(presentCols, c("PID","SID","POS","EID"))) == 2)
|| (all(is.integer(data[[presentCols[1]]])) && all(is.integer(data[[presentCols[2]]]))))) {
## additional tests above to "ensure" the two columns contain integers
keys <- .createIDs(data, cols=presentCols);
} else {
## paste the columns together
#exprStr <- paste0("paste(", paste0(paste0("data$", presentCols), collapse=", "),");");## no longer seems to work
#keys <- eval(parse(text=exprStr));
keys = apply(data[,presentCols],1,paste0,collapse=".")
}
## at this point, 'keys' is a vector
if (any(duplicated(keys))) {
#browser();return()
return(paste("The 'key' for each record is not unique.\n",
"Columns in key: ", paste(presentCols, collapse = ", "), ".\n", sep=""));
}
}
## check for numeric columns
presentCols <- intersect(numericCols, names(data));
if (length(presentCols) > 0) {
exprStr <- paste(paste("any(!is.numeric(data$", presentCols, "))", sep=""), collapse=" || ");
if (eval(parse(text=exprStr))) {
return(paste("One or more columns requires numeric values, but contains non-numerics.\n",
"Columns that must contain numerics: ", paste(presentCols, collapse = ", "),
".\n", sep=""));
}
}
## Check for increasing/descreasing POS
if (!is.null(className) && className == "PolySet") {
idx <- .createIDs(data, cols=c("PID", "SID"));
idxFirst <- which(!duplicated(idx));
idxLast <- c((idxFirst-1)[-1], length(idx));
## identify the holes
holes <- (data$POS[idxFirst] > data$POS[idxLast])
## outer/inner contour indices
idxOuter <- rep(!holes, times=((idxLast-idxFirst)+1))
idxInner <- !idxOuter;
# POS[i] < POS[i+1]?
lt <- c(data$POS[1:(nrow(data)-1)] < data$POS[2:(nrow(data))], FALSE);
## Check outer contours; change last vertex of each polygon to
## what we expect for valid outer contours
lt[idxLast] <- TRUE;
## Check for any that aren't in order
j <- any(!lt[idxOuter])
if (j) {
j <- !lt;
j[idxInner] <- FALSE;
# add 1 because it's actually the next row that break it
j <- which(j) + 1;
return(paste("POS column must contain increasing values for outer contours.\n",
"Offending rows: ", paste(j, collapse=", "), ".\n", sep=""));
}
## Check inner contours; change last vertex of each polygon to what
## we expect for valid inner contours
lt[idxLast] <- FALSE;
# check for any that aren't in order
j <- any(lt[idxInner])
if (j) {
j <- lt;
j[idxOuter] <- FALSE;
# do not add 1
j <- which(j);
return(paste("POS column must contain decreasing values for inner contours.\n",
"Offending rows: ", paste(j, collapse=", "), ".\n", sep=""));
} ## end if j
## Check for orphaned holes or hole vertices outside their solids
## Note: takes a long time if there are thousand of holes (e.g. 'rca')
## Originally reported all holes but now only reports bad hole vertices
## (any hole vertex outside the enclosing solid)
if (any(class(data) == "PolySet")) {
polys = data
polys$idx <- .createIDs(polys, cols=c("PID", "SID"));
names(holes) = unique(idx) ## assume parallelism from above
if (any(holes)) {
solids = !holes
snams = names(solids)[solids]
snums = (1:length(solids))[solids]
plist = list()
for (i in 1:length(snums)) {
ii = snums[i]
iii = snams[i]
if (iii==rev(names(solids))[1] || solids[ii+1]) next
endhole = ifelse (is.na(snums[i+1]) && rev(holes)[1], length(holes), (snums[i+1]-1))
iholes = holes[(ii+1):endhole]
if (!all(iholes)) stop ("Check '.validateData' code")
jjj = names(iholes)
jdx = idx[is.element(idx,jjj)]
isolid = polys[is.element(polys$idx,iii),]
ihole = polys[is.element(polys$idx, jjj),]
e.in.p = .is.in(ihole, isolid)
if (e.in.p$all.in) next
e.bad = e.in.p$e.out
plist[[iii]] = list()
plist[[iii]][["solid"]] = isolid
plist[[iii]][["hole.vout"]] = list()
j.bad = split(e.bad, e.bad$idx)
plist[[iii]][["hole.vout"]] = j.bad
#if (length(j.bad)>1) {browser();return()}
}
#browser();return()
if (length(plist)>0) {
assign ("solids_holes", plist, envir=.PBSmapEnv)
swo = sapply(plist, function(x){length(x$hole.vout)}) ## solids with # of orphan holes
## solids with # of holes that have vertices outside enclosing solids (originally 'solids with orphans')
swo = sapply(plist, function(x){nh=length(x$hole.vout); nv=sapply(x$hole.vout,nrow); paste0(nh, " hole",ifelse(nh>1,"s",""),": ", paste0("[",names(nv),"]=", nv,"v",collapse="; "))})
swo.txt = paste0(" Solid ",names(swo)," --> ", swo, collapse="\n")
cat(paste0("\n******* WARNING *******\n", fnam, ": Hole vertices (v) exist outside solids [PID.SID]:\n", swo.txt, "\nSee object '.PBSmapEnv$solids_holes' for details.\n\n"))
}
} ## end if any holes
} ## end if class==PolySet
#browser();return()
} ## end if className==PolySet
} else {
return(paste("The object must be either a matrix or a data frame.\n"));
}
return(data);
}
##~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~.validateData
#==============================================================================
.validateEventData <- function(EventData)
# Perform some simple tests on the object to see if it can possibly be
# an EventData object.
# If the object is invalid, returns the error message.
{
return(.validateData(EventData,
className = "EventData",
requiredCols = c("EID", "X", "Y"),
requiredAttr = NULL,
noFactorCols = c("EID", "X", "Y"),
noNACols = c("EID", "X", "Y"),
keyCols = c("EID"),
numericCols = c("EID", "X", "Y")));
}
#==============================================================================
.validateLocationSet <- function(LocationSet)
# Perform some simple tests on the object to see if it can possibly be
# an EventData object.
# If the object is invalid, returns the error message.
{
return(.validateData(LocationSet,
className = "LocationSet",
requiredCols = c("EID", "PID", "Bdry"),
requiredAttr = NULL,
noFactorCols = c("EID", "PID", "SID", "Bdry"),
noNACols = c("EID", "PID", "SID", "Bdry"),
keyCols = c("EID", "PID", "SID"),
numericCols = c("EID", "PID", "SID")));
}
#==============================================================================
.validatePolyData <- function(PolyData)
# Perform some simple tests on the object to see if it can possibly be
# a PolyData object.
# If the object is invalid, returns the error message.
{
return(.validateData(PolyData,
className = "PolyData",
requiredCols = c("PID"),
requiredAttr = NULL,
noFactorCols = c("PID", "SID"),
noNACols = c("PID", "SID"),
keyCols = c("PID", "SID"),
numericCols = c("PID", "SID")));
}
#==============================================================================
# 'polyProps' can be event properties, rather than polygon properties
.validatePolyProps <- function(polyProps, parCols = NULL)
{
if (is.null(polyProps))
return (NULL);
return (.validateData(polyProps,
className = "PolyData",
requiredCols = NULL, # PID/SID or EID
requiredAttr = NULL,
noFactorCols = parCols,
noNACols = NULL,
keyCols = NULL,
numericCols = NULL));
}
#==============================================================================
.validatePolySet <- function(polys)
# Perform some simple tests on the object to see if it can possibly be
# a PolySet object.
# If the object is invalid, returns the error message.
{
return(.validateData(polys,
className = "PolySet",
requiredCols = c("PID", "POS", "X", "Y"),
requiredAttr = NULL,
noFactorCols = c("PID", "SID", "POS", "X", "Y"),
noNACols = c("PID", "SID", "POS", "X", "Y"),
keyCols = c("PID", "SID", "POS"),
numericCols = c("PID", "SID", "POS", "X", "Y")));
}
#==============================================================================
.validateXYData <- function(xyData)
# Perform some simple tests on the object to see if it can possibly be
# an xyData object.
# If the object is invalid, returns the error message.
{
return(.validateData(xyData,
className = NULL,
requiredCols = c("X", "Y"),
requiredAttr = NULL,
noFactorCols = c("X", "Y"),
noNACols = c("X", "Y"),
keyCols = NULL,
numericCols = c("X", "Y")));
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.