R <- matrix(.5, 4, 4)
diag(R) <- 1
Npoints <- 1000
Rsq <- .40
NumDeg <- 20
V <- eigen(R)$vectors
## create b parallel to v[,3]
## rotate in the 2 - 4 plane
b <- V[,3]
bsq <- t(b) %*% R %*% b
b <- b * sqrt(Rsq/bsq)
b.i <- rcone(R, Rsq,b, V[,2], V[,4], deg = NumDeg, Npoints)
t(b.i[,1]) %*% R %*% b.i[,1]
t(b.i[,25]) %*% R %*% b.i[,25]
Run the code above in your browser using DataLab