[R] rgl: cylinder3d() with elliptical cross-section

[Duncan Murdoch]

On 10/03/2012 9:44 AM, Duncan Murdoch wrote:
> My first reply to this went privately, by accident.  I've done a little
> editing to it, but mainly this is for the archives.
>
> On 12-03-09 2:36 PM, Michael Friendly wrote:
> >  For a paper dealing with generalized ellipsoids, I want to illustrate in
> >  3D an ellipsoid that is unbounded
> >  in one dimension, having the shape of an infinite cylinder along, say,
> >  z, but whose cross-section in (x,y)
> >  is an ellipse, say, given by the 2x2 matrix cov(x,y).
> >
> >  I've looked at rgl:::cylinder3d, but don't see any way to make it
> >  accomplish this.  Does anyone have
> >  any ideas?
>
> rgl has no way to display curved surfaces that are unbounded.  (It has
> lines and planes that adapt to the viewport.)  So you would need to make
> a finite cylinder, and it will be up to you to choose how big to make it.
>
> The cylinder3d() function can do that, but it's not very good at
> cylinders that are straight. (This is a little embarrassing...) It sets
> up a local coordinate system based on the curvature, but if there are no
> curves, it fails, and you have to supply your own coordinates.

This got forgotten for a while, but is now fixed on R-forge.  If you get 
rev 881 or newer of rgl, it can draw straight cylinders.

Duncan Murdoch

>
> So here's how I would do what you want:
>
> center<- cbind(0, 0, 1:10)  # cylinder centered on points (0, 0, z)
> e2<- cbind(1, 0, rep(0, 10)) # define the normal vectors
> cyl<- cylinder3d(center, e2=e2)
>
> # Now you have an octagonal cylinder.  Use the sides arg to cylinder3d
> # if it doesn't end up smooth enough, but in most cases I've seen, 8
> # is sufficient.
>
> # Define a transformation to the x and y coordinates to give the
> # elliptical shape; use it as the
> # top left 2x2 matrix of a 3x3 matrix
> xfrm<- matrix( c(2, 1, 0,
>                                  1, 3, 0,
>                                  0, 0, 1), 3,3, byrow=TRUE)
> cyl<- transform3d(cyl, xfrm)
> cyl<- addNormals(cyl)  # this makes it shade smoothly
> shade3d(cyl, col="green")
> decorate3d()  # show some axes for scale
>
> Duncan Murdoch