# [R] orthographic projection of ellipsoids

baptiste auguie ba208 at exeter.ac.uk
Sun Oct 26 14:26:11 CET 2008

```Dear list,

I've generated a list of 3D coordinates representing ellipsoids in
arbitrary orientations. I'm now trying to obtain a 2D projection of
the scene, that is to draw the silhouette of each object on a plane
(x,y). The only way I could think of is to compute the convex hull of
the (x,y) coordinates of each object and use this as the outline of
the object. This is clearly not very efficient or satisfying.

I think I'm on the wrong track from the start. Is there an obvious
analytical parametrisation of such projections? Any comments are
welcome.

Many thanks,

baptiste

>
> rotM3d <- function(theta=0, phi=0, psi=0){ # 3D rotation matrix
> 	a11 <- cos(psi)*cos(phi) - cos(theta)*sin(phi)*sin(psi)
> 	a12 <- cos(psi)*sin(phi) + cos(theta)*cos(phi)*sin(psi)
> 	a13 <- sin(psi)*sin(theta)
> 	a21 <- -sin(psi)*cos(phi) - cos(theta)*sin(phi)*cos(psi)
> 	a22 <- -sin(psi)*sin(phi) + cos(theta)*cos(phi)*cos(psi)
> 	a23 <- cos(psi)*sin(theta)
> 	a31 <- sin(theta)*sin(phi)
> 	a32 <- -sin(theta)*cos(phi)
> 	a33 <- cos(theta)
> 	matrix(c(a11, a12, a13, a21, a22, a23, a31, a32, a33), ncol=3)
> }
> rotM3d() # I
>
> ellipsoid <-  # idea borrowed from a post in the R-mailing list
> (John Fox i think)
> function(x=0, y=0, z=0, radius=1, shape=diag(c(10, 2, 2)),theta=0,
> phi=0, psi=0,  segments=11) {
> 			 angles <- (0:segments)*2*pi/segments
> 			ecoord2 <- function(p) {
> 			    c(cos(p)*sin(p), sin(p)*sin(p), cos(p))
> 			}
> unit.sphere <- t(apply(expand.grid(angles, angles), 1, ecoord2))
> xyz <- t(c(x, y, z) + radius * rotM3d(theta, phi, psi)%*
> %t(unit.sphere %*% chol(shape)))
> chull(x=xyz[, 1], y=xyz[, 2])->points
> mdf <- data.frame(x=xyz[points, 1], y=xyz[points, 2])
> polygon(mdf, col=hcl(h = 0, c = 35, l = 85, 0.5))
> invisible(xyz)
> }
>
>
> xx <- seq(-5, 5, len=10)
> xy <- expand.grid(xx, xx)
>
> xy.jit <- apply(xy, 2, jitter, amount=0.4)
>
> par(mar=c(0, 0, 0, 0))
> plot(xy.jit, t="n", axes=F, xlab="", ylab="")
>
> x <- xy.jit[, 1]
> y <- xy.jit[, 2]
>
> twist <- pi*y/max(abs(y)) * rep(1, length(y))
> tilt <- pi*x/max(abs(x)) * rep(1, length(x))
> b.quiet <- mapply(ellipsoid,
> 	theta=twist, psi=tilt,x=x, y=y,  SIMPLIFY=F, radius=0.15)

_____________________________

Baptiste Auguié

School of Physics
University of Exeter