Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

How do I draw a 3D arrow in R?

I am using the rgl package to create 3D plots of my data. For some reasons (3D PCA biplots) I need vectors -- a line segment with an arrow. And I'm stuck, because I want to have 3D cones as arrow heads.

Somehow, I cannot wrap my senile mind around the geometry of the problem. Say, I would draw the vector with

segments3d( rbind( c( 0, 0, 0 ), c( 3, 3, 3 ) ) )

that is, a vector from the origin of the user coordinate system to [3,3,3].

I would like to create a cone with the tip at [3,3,3]. The base of the cone can be formed with a circle. Drawing a circle on the xz plane (perpendicular to the y plane) with radius r is easy:

n <- 10
sin.t <- sin( seq( 0, 2 * pi, len= n ) )
cos.t <- cos( seq( 0, 2 * pi, len= n ) )
r <- 0.1 
xv <- x + r * sin.t
yv <- rep( y, n )
zv <- z + r * cos.t

but how do I now transform these points such that the circle is now perpendicular to the vector? And its center 0.2 from the tip along the vectors direction? Once I have this transformation, I will draw triangles with the triangles3d function, each triangle having one corner at the tip and two vertices within the points of the circle.

This is basic maths, and I know the 18 year old me would not have a problem (or even a 28 year old me). Any hook (as opposed to fish) would be appreciated.

like image 732
January Avatar asked Mar 01 '13 20:03

January


People also ask

How do you make an arrow in R?

To create an arrow R, we can use plot function and arrows function. We just need to understand all the coordinate values that should be passed inside the arrows function. For example, if we have two vectors that contains values from 1 to 10 then the arrow can be created by using arrows function as arrows(1,1,10,10).


1 Answers

In the demos for rgl there is a cone3d function. It takes the base and the tip seperaately. In any event you could do something like this:

vec=rbind( c( 0, 0, 0 ), c( 3, 3, 3 ) )
segments3d( vec )


cone3d(base=vec[2,]-(vec[1,]+vec[2,]/6), 
     #this makes the head go 1/6th the length of the arrow
       rad=0.5,
       tip=vec[2,],
       col="blue",
       front="lines",
       back="lines")

Here is the cone3d function:

   cone3d <- function(base=c(0,0,0),tip=c(0,0,1),rad=1,n=30,draw.base=TRUE,qmesh=FALSE,
                    trans = par3d("userMatrix"), ...) {
   ax <- tip-base
   if (missing(trans) && !rgl.cur()) trans <- diag(4)
   ### is there a better way?
   if (ax[1]!=0) {
     p1 <- c(-ax[2]/ax[1],1,0)
     p1 <- p1/sqrt(sum(p1^2))
     if (p1[1]!=0) {
       p2 <- c(-p1[2]/p1[1],1,0)
       p2[3] <- -sum(p2*ax)
       p2 <- p2/sqrt(sum(p2^2))
     } else {
       p2 <- c(0,0,1)
     }
   } else if (ax[2]!=0) {
     p1 <- c(0,-ax[3]/ax[2],1)
     p1 <- p1/sqrt(sum(p1^2))
     if (p1[1]!=0) {
       p2 <- c(0,-p1[3]/p1[2],1)
       p2[3] <- -sum(p2*ax)
       p2 <- p2/sqrt(sum(p2^2))
     } else {
       p2 <- c(1,0,0)
     }
   } else {
     p1 <- c(0,1,0); p2 <- c(1,0,0)
   }
   degvec <- seq(0,2*pi,length=n+1)[-1]
   ecoord2 <- function(theta) {
     base+rad*(cos(theta)*p1+sin(theta)*p2)
   }
   i <- rbind(1:n,c(2:n,1),rep(n+1,n))
   v <- cbind(sapply(degvec,ecoord2),tip)
   if (qmesh) 
     ## minor kluge for quads -- draw tip twice
     i <- rbind(i,rep(n+1,n))
   if (draw.base) {
     v <- cbind(v,base)
     i.x <- rbind(c(2:n,1),1:n,rep(n+2,n))
     if (qmesh)  ## add base twice
       i.x <-  rbind(i.x,rep(n+2,n))
     i <- cbind(i,i.x)
   }
   if (qmesh) v <- rbind(v,rep(1,ncol(v))) ## homogeneous
   if (!qmesh)
     triangles3d(v[1,i],v[2,i],v[3,i],...)
   else
     return(rotate3d(qmesh3d(v,i,material=...), matrix=trans))
 }     
like image 146
Seth Avatar answered Oct 17 '22 13:10

Seth