[R] plotting an isosurface on a 3d plot

[Duncan Murdoch]

On 13/06/2019 4:32 p.m., Duncan Murdoch wrote:
> On 13/06/2019 12:47 p.m., ravi via R-help wrote:
>> Hi,I want to plot a surface joining a circle on a plane with another circle (a little offset w.r.t. the first one) on a parallel plane. This is a very simplified version of my problem.
>> I will explain with the following code :
>> # First, I plot the two circlesorig1 <- 0.4;radius1=0.3;theta <- seq(0,2*pi,pi/8)x1 <- orig1+radius1*cos(theta)y1 <- orig1+radius1*sin(theta)
>> orig2 <- 0.7;radius2=0.2;theta <- seq(0,2*pi,pi/8)x2 <- orig2+radius2*cos(theta)y2 <- orig2+radius2*sin(theta)
>> plot(x1,y1,type='b',col="red",xlim=c(0,1),ylim=c(0,1),asp=1)lines(x2,y2,type="b",col="blue")
>> #### the z coordinates on two parallel plnesz1 <- rep(0,length(x1))z2 <- rep(1,length(x2))
>>
>> I have tried to plot my 3d figure using the packages plot3D, misc3d, rgl etc. All these packages require z as a matrix. In my case, all of x,y and z are vectors. How can I plot this surface? In addition, I would like to add similiar surfaces to the same plot.
> 
> It's not true that rgl requires z as a matrix:  that's one possibility,
> but there are others.  Here's one way to do what you want.
> 
> # First, compute the two circles
> orig1 <- 0.4;radius1=0.3;theta <- seq(0,2*pi,pi/8)
> x1 <- orig1+radius1*cos(theta)
> y1 <- orig1+radius1*sin(theta)
> orig2 <- 0.7;radius2=0.2
> x2 <- orig2+radius2*cos(theta)
> y2 <- orig2+radius2*sin(theta)
> 
> #### the z coordinates on two parallel plnes
> z1 <- rep(0,length(x1))
> z2 <- rep(1,length(x2))
> 
> library(rgl)
> p1 <- cbind(x1,y1,z1)
> plot3d(p1, col="red")
> p2 <- cbind(x2, y2, z2)
> points3d(p2, col="blue")
> 
> quads <- matrix(numeric(), 0, 3)
> for (i in seq_along(x1)[-1]) {
>     quads <- rbind(quads, p1[i-1,], p1[i,], p2[i,], p2[i-1,])
> }
> quads3d(quads, col = "green")
> 
> There are more efficient ways to do this (the for loop would be slow if
> you had bigger vectors), but this should give you the idea.

Here's one of the more efficient ways:

cylinder3d(rbind(c(0.4, 0.4, 0), c(0.7, 0.7, 1)),
            radius = c(0.3, 0.2),
            e1 = rbind(c(0,0,1), c(0,0,1)), sides=20) %>%
            addNormals() %>%
            shade3d(col = "green")

Duncan Murdoch