Trisurf Plots in R using Plotly

June 1, 2016
By

(This article was first published on R Modern Data, and kindly contributed to R-bloggers)

In this post well show how to create Triangular Surface Plots in R. This post is based on timelyportfolios gist.

Moebius Strip

library(plotly)
library(geometry)

g <- expand.grid(
  u = seq(0, 2 * pi, length.out = 24),
  v = seq(-1, 1, length.out = 8)
)
tp <- 1 + 0.5 * g$v * cos(g$u / 2)
m <- matrix(
  c(tp * cos(g$u), tp * sin(g$u), 0.5 * g$v * sin(g$u / 2)), 
  ncol = 3, dimnames = list(NULL, c("x", "y", "z"))
)

# the key though is running delaunayn on g rather than m
d <- delaunayn(g)
td <- t(d)
#  but using m for plotting rather than the 2d g

# define layout options
axs <- list(
  backgroundcolor="rgb(230,230,230)",
  gridcolor="rgb(255,255,255)",
  showbackground=TRUE,
  zerolinecolor="rgb(255,255,255"
)

# now figure out the colormap
#  start by determining the mean of z for each row
#  of the Delaunay vertices
zmean <- apply(d, MARGIN=1, function(row){mean(m[row,3])})

library(scales)
# result will be slighlty different
#  since colour_ramp uses CIELAB instead of RGB
#  could use colorRamp for exact replication
facecolor = colour_ramp(
  brewer_pal(palette="RdBu")(9)
)(rescale(x=zmean))


plot_ly(
  x = m[, 1], y = m[, 2], z = m[, 3],
  # JavaScript is 0 based index so subtract 1
  i = d[, 1]-1, j = d[, 2]-1, k = d[, 3]-1,
  facecolor = facecolor,
  type = "mesh3d"
) %>%
  layout(
    title="Moebius band triangulation",
    scene=list(xaxis=axs,yaxis=axs,zaxis=axs),
    aspectratio=list(x=1,y=1,z=0.5)
  )

2D Surface over a disk

n <- 12
h <- 1/(n-1)
r = seq(h, 1, length.out=n)
theta = seq(0, 2*pi, length.out=36)

g <- expand.grid(r=r, theta=theta)

x <- c(g$r * cos(g$theta),0)
y <- c(g$r * sin(g$theta),0)
z <- sin(x*y)

m <- matrix(
  c(x,y,z), 
  ncol = 3,
  dimnames = list(NULL, c("x", "y", "z"))
)

tri <- delaunayn(m[,1:2])

# now figure out the colormap
zmean <- apply(tri,MARGIN=1,function(row){mean(m[row,3])})

library(scales)
library(rje)
facecolor = colour_ramp(
  cubeHelix(12)
)(rescale(x=zmean))

plot_ly(
  x=x, y=y, z=z,
  i=tri[,1]-1, j=tri[,2]-1, k=tri[,3]-1,
  facecolor=facecolor,
  type="mesh3d"
) %>%
  layout(
    title="Triangulated surface",
    scene=list(
      xaxis=axs,
      yaxis=axs,
      zaxis=axs,
      camera=list(
        eye=list(x=1.75,y=-0.7,z=0.75)
      )
    ),
    aspectratio=list(x=1,y=1,z=0.5)
  )

Chopper from python

library(geomorph)
plyFile <- 'http://people.sc.fsu.edu/~jburkardt/data/ply/chopper.ply'
dest <- basename(plyFile)
if (!file.exists(dest)) {
  download.file(plyFile, dest)
}

mesh <- read.ply(dest)
# see getS3method("shade3d", "mesh3d") for details on how to plot 

# plot point cloud
x <- mesh$vb["xpts",]
y <- mesh$vb["ypts",]
z <- mesh$vb["zpts",]
m <- matrix(c(x,y,z), ncol=3, dimnames=list(NULL,c("x","y","z")))

# now figure out the colormap
zmean <- apply(t(mesh$it),MARGIN=1,function(row){mean(m[row,3])})

library(scales)
facecolor = colour_ramp(
  brewer_pal(palette="RdBu")(9)
)(rescale(x=zmean))

plot_ly(
  x = x, y = y, z = z,
  i = mesh$it[1,]-1, j = mesh$it[2,]-1, k = mesh$it[3,]-1,
  facecolor = facecolor,
  type = "mesh3d"
)

To leave a comment for the author, please follow the link and comment on their blog: R Modern Data.

R-bloggers.com offers daily e-mail updates about R news and tutorials on topics such as: Data science, Big Data, R jobs, visualization (ggplot2, Boxplots, maps, animation), programming (RStudio, Sweave, LaTeX, SQL, Eclipse, git, hadoop, Web Scraping) statistics (regression, PCA, time series, trading) and more...



If you got this far, why not subscribe for updates from the site? Choose your flavor: e-mail, twitter, RSS, or facebook...

Comments are closed.

Search R-bloggers


Sponsors

Never miss an update!
Subscribe to R-bloggers to receive
e-mails with the latest R posts.
(You will not see this message again.)

Click here to close (This popup will not appear again)