Dear list,
I'm trying to visualise some ellipsoidal shapes in 3D. Their position,
axes, and angular orientation can be arbitrary. I saw an ellipse3d
function in rgl; however it is heavily oriented towards the
statistical concept of ellipse of confidence, whilst I am just
concerned with the geometrical object.
Below is my current implementation. It is quite slow when creating
many shapes, so I'm asking for suggestions (I'm new to rgl),
## helping function, probably already in rgl but nevermind
euler <- function(alpha=0, beta=0, gamma=0){
Ra <- matrix(c(cos(alpha), sin(alpha), 0, -sin(alpha), cos(alpha),
0, 0, 0, 1), ncol=3, byrow=T)
Rb <- matrix(c(1, 0, 0, 0, cos(beta), sin(beta), 0, -sin(beta),
cos(beta)), ncol=3, byrow=T)
Rc <- matrix(c(cos(gamma), sin(gamma), 0, -sin(gamma), cos(gamma),
0, 0, 0, 1), ncol=3, byrow=T)
return(Rc%*%Rb%*%Ra)
}
require(rgl)
## create one ellipsoid at a given location, with axes a,b,c and Euler
angles phi, theta, psi
rgl.ellipsoid <- function (x=0,y=0,z=0, a = 1,b=1,c=1, phi=0,theta=0,psi=0,
subdivide = 3, smooth = TRUE, ...)
{
sphere <- subdivision3d(cube3d(...), subdivide)
class(sphere) <- c("mesh3d","shape3d")
norm <- sqrt(sphere$vb[1, ]^2 + sphere$vb[2, ]^2 + sphere$vb[3,
]^2)
for (i in 1:3) sphere$vb[i, ] <- sphere$vb[i, ]/norm
sphere$vb[4, ] <- 1
sphere$normals <- sphere$vb
result <- scale3d(sphere, a,b,c)
rotM <- euler(phi,theta,psi)
result <- rotate3d(result,matrix=rotM)
result <- translate3d(result, x,y,z)
invisible(result)
}
## loop over the specification of a cluster
rgl.ellipsoids <- function(positions, sizes, angles,...){
N <- NROW(positions)
ll <- lapply(seq(1,N), function(ii)
rgl.ellipsoid(positions[ii,1],positions[ii,2],positions[ii,3],
sizes[ii,1],sizes[ii,2],sizes[ii,3],
angles[ii,1],angles[ii,2],angles[ii,3], ...))
shapelist3d(ll,...)
}
N <- 20
set.seed(123)
positions <- matrix(rnorm(3*N), ncol=3)
sizes <- matrix(runif(3*N, 0.01, 0.05), ncol=3)
angles <- matrix(runif(3*N, 0, 2*pi), ncol=3)
rgl.ellipsoids(positions, sizes, angles, col= 1:N)
Best regards,
baptiste
sessionInfo()
R version 2.12.1 Patched (2010-12-30 r53895)
Platform: i386-apple-darwin9.8.0 (32-bit)
locale:
[1] en_GB.UTF-8/en_GB.UTF-8/C/C/en_GB.UTF-8/en_GB.UTF-8
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] rgl_0.92.798
loaded via a namespace (and not attached):
[1] tools_2.12.1
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
