# R/rcsr.R In HiveR: 2D and 3D Hive Plots for R

#### Documented in rcsr

```#' Compute the Details of a 3D Spline for a Hive Plot Edge
#'
#' This is a wild bit of trigonometry!  Three points in 3D space, two ends and
#' an control point, are rotated into 2D space.  Then a spline curve is
#' computed.  This is necessary because spline curves are only defined in
#' \code{R} as 2D objects.  The new collection of points, which is the complete
#' spline curve and when drawn will be the edge of a hive plot, is rotated back
#' into the original 3D space. \code{rcsr} stands for rotate, compute spline,
#' rotate back.
#'
#' See the code for exactly how the function works.  Based upon the process
#' described at \url{http://www.fundza.com/mel/axis_to_vector/index.html}
#' Timing tests show this function is fast and scales linearly (i.e. 10x more
#' splines to draw takes 10x more time).  Roughly 3 seconds were required to
#' draw 1,000 spline curves in my testing.
#'
#' @param p0 A triple representing one end of the final curve (x, y, z).
#'
#' @param cp A triple representing the control point used to compute the final
#' curve (x, y, z).
#'
#' @param p1 A triple representing the other end of the final curve (x, y, z).
#'
#' @return A 3 column matrix with the x, y and z coordinates to be plotted to
#' create a hive plot edge.
#'
#' @author Bryan A. Hanson, DePauw University. \email{[email protected]@depauw.edu}
#'
#' @keywords utilities 3D spline
#'
#' @export rcsr
#'
#' @importFrom stats spline
#'
#' @examples
#'
#' # This is a lengthy example to prove it works.
#' # Read it and then copy the whole thing to a blank script.
#' # Parts of it require rgl and are interactive.
#' # So none of the below is run during package build/check.
#'
#' ### First, a helper function
#' \dontrun{
#'
#' drawUnitCoord <- function() {
#'
#' 	# Simple function to draw a unit 3D coordinate system
#'
#' 	# Draw a Coordinate System
#'
#' 	r <- c(0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1) # in polar coordinates
#' 	theta <- c(0, 0, 0, 90, 0, 180, 0, 270, 0, 0, 0, 0)  # start, end, start, end
#' 	phi <- c(0, 90, 0, 90, 0, 90, 0, 90, 0, 0, 0, 180)
#' 	cs <- data.frame(radius = r, theta, phi)
#' 	ax.coord <- sph2cart(cs)
#'
#' 	segments3d(ax.coord, col = "gray", line_antialias = TRUE)
#' 	points3d(x = 0, y = 0, z = 0, color = "black", size = 4,
#' 		point_antialias = TRUE) # plot origin
#'
#' 	# Label the axes
#'
#' 	r <- c(1.1, 1.1, 1.1, 1.1, 1.1, 1.1) # in polar coordinates
#' 	theta <- c(0, 90, 180, 270, 0, 0)
#' 	phi <- c(90, 90, 90, 90, 0, 180)
#' 	l <- data.frame(radius = r, theta, phi)
#' 	lab.coord <- sph2cart(l)
#' 	text3d(lab.coord, texts = c("+x", "+y", "-x", "-y", "+z", "-z"))
#'
#' 	}
#'
#' ###  Now, draw a reference coordinate system and demo the function in it.
#'
#' drawUnitCoord()
#'
#' ### Draw a bounding box
#'
#' box <- data.frame(
#' 	x = c(1, -1, 1, 1, 1, 1, 1, 1, 1, -1, -1, -1, 1, 1, 1, -1, 1, -1, -1, -1, -1, -1, -1, -1),
#' 	y = c(1, 1, 1, 1, 1, -1, 1, -1, 1, 1, 1, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, 1, -1, 1),
#' 	z = c(1, 1, 1, -1, 1, 1, -1, -1, -1, -1, 1, -1, 1, -1, 1, 1, -1, -1, -1, 1, 1, 1, -1, -1))
#'
#' segments3d(box\$x, box\$y, box\$z, line_antialias = TRUE, col = "red")
#'
#' ### Draw the midlines defining planes
#'
#' mid <- data.frame(
#' 	x = c(0, 0, 0, 0, 0, 0, 0, 0, 1, -1, -1, -1, -1, 1, 1, 1, 1, -1, -1, -1, -1, 1, 1, 1),
#' 	y = c(-1, -1, -1, 1, 1, 1, 1, -1, 0, 0, 0, 0, 0, 0, 0, 0, -1, -1, -1, 1, 1, 1, 1, -1),
#' 	z = c(-1, 1, 1, 1, 1, -1, -1, -1, -1, -1, -1, 1, 1, 1, 1, -1, 0, 0, 0, 0, 0, 0, 0, 0))
#'
#' segments3d(mid\$x, mid\$y, mid\$z, line_antialias = TRUE, col = "blue")
#'
#' ### Generate two random points
#'
#' p <- runif(6, -1, 1)
#'
#' # Special case where p1 is on z axis
#' # Uncomment line below to demo
#' #p[4:5] <- 0
#'
#' p0 <- c(p[1], p[2], p[3])
#' p1 <- c(p[4], p[5], p[6])
#'
#' ### Draw the pts, label them, draw vectors to those pts from origin
#'
#' segments3d(x = c(0, p[1], 0, p[4]),
#' 	y = c(0, p[2], 0, p[5]),
#' 	z = c(0, p[3], 0, p[6]),
#' 	line_antialias = TRUE, col = "black", lwd = 3)
#'
#' points3d(x = c(p[1], p[4]),
#' 	y = c(p[2], p[5]),
#' 	z = c(p[3], p[6]),
#' 	point_antialias = TRUE, col = "black", size = 8)
#'
#' text3d(x = c(p[1], p[4]),
#' 	y = c(p[2], p[5]),
#' 	z = c(p[3], p[6]),
#' 	col = "black", texts = c("p0", "p1"), adj = c(1,1))
#'
#' ### Locate control point
#' ### Compute and draw net vector from origin thru cp
#' ### Connect p0 and p1
#'
#' s <- p0 + p1
#' segments3d(x = c(0, s[1]), y = c(0, s[2]), z = c(0, s[3]),
#' 	line_antialias = TRUE, col = "grey", lwd = 3)
#'
#' segments3d(x = c(p[1], p[4]), # connect p0 & p1
#' 	y = c(p[2], p[5]),
#' 	z = c(p[3], p[6]),
#' 	line_antialias = TRUE, col = "grey", lwd = 3)
#'
#' cp <- 0.6*s # Now for the control point
#'
#' points3d(x = cp[1], # Plot the control point
#' 	y = cp[2],
#' 	z = cp[3],
#' 	point_antialias = TRUE, col = "black", size = 8)
#'
#' text3d(x = cp[1], # Label the control point
#' 	y = cp[2],
#' 	z = cp[3],
#' 	texts = c("cp"), col = "black", adj = c(1,1))
#'
#' ### Now ready to work on the spline curve
#'
#' n2 <- rcsr(p0, cp, p1) # Compute the spline
#'
#' lines3d(x = n2[,1], y = n2[,2], z = n2[,3],
#' 	line_antialias = TRUE, col = "blue", lwd = 3)
#'
#' ### Ta-Da!!!!!
#' }
#'
rcsr <- function(p0, cp, p1) {

# Bryan Hanson, DePauw Univ, April 2011

# Function to take 3 points in 3d and compute a spline curve
# using them.  The strategy is to rotate into a 2d system,
# figure the spline, and then rotate back to the orig 3d system.

# p0 and p1 are end points, cp is the control point

# rcsr = rotate, compute spline, rotate back

# Local copies of functions borrowed from RFOC

rotx3 <- function (deg) {
r = diag(3)
r[3, 3] = r[2, 2]
r[3, 2] = -r[2, 3]
return(r)
}

roty3 <- function (deg) {
r = diag(3)
r[3, 3] = r[1, 1]
r[1, 3] = -r[3, 1]
return(r)
}

rotz3 <- function (deg) {
r = diag(3)
r[2, 2] = r[1, 1]
r[2, 1] = -r[1, 2]
return(r)
}

m <- matrix(c(0, 0, 0, p0, cp, p1), nrow = 4, byrow = TRUE)

# Align p0 with the +y axis by rotating around x and z axes
# see www.fundza.com/mel/axis_to_vector/align_axis_to_vector.html
# Other points follow using rotation matrix

xyL <- sqrt(m[2,1]^2 + m[2,2]^2)
vL <- sqrt(sum((m[2,]^2)))

if (xyL == 0) { # have to catch when p1 is on z axis
zA <- ifelse(m[2,1] > 0, 90, -90)
} else {
zA <-acos(m[2,2]/xyL)
}
zA <- zA*180/pi
zA <- rotz3(zA)

xA <- acos(xyL/vL)
xA <- xA*180/pi
xA <- rotx3(xA)

n1 <- t(zA %*% t(m)) # order of rotations matter!
n2 <- t(xA %*% t(n1))

# Now rotate p1 to +x+y plane by rotating around y axis

yzL <- sqrt(n2[4,1]^2 + n2[4,3]^2)

yA <-acos(n2[4,3]/yzL)
yA <- 90 + (yA*180/pi)
yA <- roty3(yA)

n3 <- t(yA %*% t(n2))

# Compute spline curve
n3 <- n3[-1,] # remove origin point
sp <- stats::spline(n3[,1], n3[,2], n = 25)
sp <- matrix(c(sp\$x, sp\$y, rep(0, length(sp\$x))), ncol = 3) # add back the z coord (all 0's)

# Now reverse the transformations back to the original 3d space

sp <- t(t(yA) %*% t(sp))
sp <- t(t(xA) %*% t(sp))
sp <- t(t(zA) %*% t(sp))
sp
}
```

## Try the HiveR package in your browser

Any scripts or data that you put into this service are public.

HiveR documentation built on May 2, 2019, 2:08 a.m.