Here is an approach using kernel density estimation and the contour3d
function from misc3d
. I played around until I found a value for levels
that worked decently. It's not perfectly precise, but you may be able to tweak things to get a better, more accurate surface. If you have more than 8GB of memory, then you may be able to increase n
beyond what I did here.
library(rgl)
library(misc3d)
library(onion); data(bunny)
# the larger the n, the longer it takes, the more RAM you need
bunny.dens <- kde3d(bunny[,1],bunny[,2],bunny[,3], n=150,
lims=c(-.1,.2,-.1,.2,-.1,.2)) # I chose lim values manually
contour3d(bunny.dens$d, level = 600,
color = "pink", color2 = "green", smooth=500)
rgl.viewpoint(zoom=.75)
The image on the right is from the bottom, just to show another view.
You can use a larger value for n
in kde3d
but it will take longer, and you may run out of RAM if the array becomes too large. You could also try a different bandwidth (default used here). I took this approach from Computing and Displaying Isosurfaces in R - Feng & Tierney 2008.
Very similar isosurface approach using the Rvcg
package:
library(Rvcg)
library(rgl)
library(misc3d)
library(onion); data(bunny)
bunny.dens <- kde3d(bunny[,1],bunny[,2],bunny[,3], n=150,
lims=c(-.1,.2,-.1,.2,-.1,.2)) # I chose lim values manually
bunny.mesh <- vcgIsosurface(bunny.dens$d, threshold=600)
shade3d(vcgSmooth(bunny.mesh,"HC",iteration=3),col="pink") # do a little smoothing
Since it's a density estimation based approach, we can get a little more out of it by increasing the density of the bunny. I also use n=400
here. The cost is a significant increase in computation time, but the resulting surface is a hare better:
bunny.dens <- kde3d(rep(bunny[,1], 10), # increase density.
rep(bunny[,2], 10),
rep(bunny[,3], 10), n=400,
lims=c(-.1,.2,-.1,.2,-.1,.2))
bunny.mesh <- vcgIsosurface(bunny.dens$d, threshold=600)
shade3d(vcgSmooth(bunny.mesh,"HC",iteration=1), col="pink")
Better, more efficient surface reconstruction methods exist (e.g. power crust, Poisson surface reconstruction, ball-pivot algorithm), but I don't know that any have been implemented in R, yet.
Here's a relevant Stack Overflow post with some great information and links to check out (including links to code): robust algorithm for surface reconstruction from 3D point cloud?.