R + plotly：革命实体

Question

我有一个函数r(x)，我想围绕x 轴旋转以获得一个solid of revolution，我想使用add_surface 将其添加到现有的plot_ly 图中（由x 着色） .

这是一个例子：

library(dplyr)
library(plotly)

# radius depends on x
r <- function(x) x^2

# interval of interest
int <- c(1, 3)

# number of points along the x-axis
nx <- 20

# number of points along the rotation
ntheta <- 36

# set x points and get corresponding radii
coords <- data_frame(x = seq(int[1], int[2], length.out = nx), r = r(x))

# for each x: rotate r to get y and z coordinates
# edit: ensure 0 and pi are both amongst the angles used
coords %<>%
  rowwise() %>%
  do(data_frame(x = .$x, r = .$r,
                theta = seq(0, pi, length.out = ntheta / 2 + 1) %>%
                c(pi + .[-c(1, length(.))]))) %>%

  ungroup %>%
  mutate(y = r * cos(theta), z = r * sin(theta))

# plot points to make sure the coordinates define the desired shape
coords %>%
  plot_ly(x = ~x, y = ~y, z = ~z, color = ~x) %>%
  add_markers()

如何将上述点指示的形状生成为plotly 表面（理想情况下两端都是开放的）？

---

编辑（1）：

这是我迄今为止最好的尝试：

# get all x & y values used (sort to connect halves on the side)
xs <-
  unique(coords$x) %>%
  sort
ys <-
  unique(coords$y) %>%
  sort

# for each possible x/y pair: get z^2 value
coords <-
  expand.grid(x = xs, y = ys) %>%
  as_data_frame %>%
  mutate(r = r(x), z2 = r^2 - y^2)

# format z coordinates above x/y plane as matrix where columns
# represent x and rows y
zs <- matrix(sqrt(coords$z2), ncol = length(xs), byrow = TRUE)

# format x coordiantes as matrix as above (for color gradient)
gradient <-
  rep(xs, length(ys)) %>%
  matrix(ncol = length(xs), byrow = TRUE)
  
# plot upper half of shape as surface
p <- plot_ly(x = xs, y = ys, z = zs, surfacecolor = gradient,
             type = "surface", colorbar = list(title = 'x'))

# plot lower have of shape as second surface
p %>%
  add_surface(z = -zs, showscale = FALSE)

虽然这给出了所需的形状，

1. 它在x/y 平面附近有“锋利的牙齿”。
2. 半部分不接触。（通过在theta 向量中包含0 和pi 解决）
3. 我不知道如何用x 代替z 给它上色（尽管到目前为止我并没有对此进行过多研究）。（由gradient 矩阵解决)

---

编辑（2）：

这是使用单个表面的尝试：

# close circle in y-direction
ys <- c(ys, rev(ys), ys[1])

# get corresponding z-values
zs <- rbind(zs, -zs[nrow(zs):1, ], zs[1, ])

# as above, but for color gradient
gradient <-
  rbind(gradient, gradient[nrow(gradient):1, ], gradient[1, ])

# plot single surface
plot_ly(x = xs, y = ys, z = zs, surfacecolor = gradient,
        type = "surface", colorbar = list(title = 'x'))

令人惊讶的是，虽然这应该连接与x / y 平面正交的两半以创建完整的形状， 它仍然遭受与上述解决方案相同的“剃刀齿”效果：

---

编辑（3）：

原来缺少部分是由于z-values 在接近 0 时为 NaN：

# color points 'outside' the solid purple
gradient[is.nan(zs)] <- -1

# show those previously hidden points
zs[is.nan(zs)] <- 0

# plot exactly as before
plot_ly(x = xs, y = ys, z = zs, surfacecolor = gradient,
        type = "surface", colorbar = list(title = 'x'))

这可能是由于 r^2 和 y 太接近时减法的数值不稳定导致的，导致 sqrt 的负输入，而实际输入仍然是非负的。

这与数值问题无关，因为即使考虑到 +-4“接近”为零，也无法完全避免“剃刀齿”效应：

# re-calculate z-values rounding to zero if 'close'
eps <- 4
zs <- with(coords, ifelse(abs(z2) < eps, 0, sqrt(z2))) %>%
      matrix(ncol = length(xs), byrow = TRUE) %>%
      rbind(-.[nrow(.):1, ], .[1, ])

# plot exactly as before
plot_ly(x = xs, y = ys, z = zs, surfacecolor = gradient,
        type = "surface", colorbar = list(title = 'x'))

Answer 1

有趣的问题，我一直在努力使用表面密度来改进您的解决方案。你可以通过分层多行来做一个黑客，这对这个很好，例如仅对原始数据进行更改，例如使用更多的 x 点：nx 到 1000，并将 add_markers 更改为 add_lines。可能无法扩展，但适用于这种大小的数据:)

library(dplyr)
library(plotly)

# radius depends on x
r <- function(x) x^2

# interval of interest
int <- c(1, 3)

# number of points along the x-axis
nx <- 1000

# number of points along the rotation
ntheta <- 36

# set x points and get corresponding radii
coords <- data_frame(x = seq(int[1], int[2], length.out = nx), r = r(x))

# for each x: rotate r to get y and z coordinates
# edit: ensure 0 and pi are both amongst the angles used
coords %<>%
  rowwise() %>%
  do(data_frame(x = .$x, r = .$r,
                theta = seq(0, pi, length.out = ntheta / 2 + 1) %>%
                  c(pi + .[-c(1, length(.))]))) %>%

  ungroup %>%
  mutate(y = r * cos(theta), z = r * sin(theta))

# plot points to make sure the coordinates define the desired shape
coords %>%
  plot_ly(x = ~x, y = ~y, z = ~z, color = ~x) %>%
  add_lines()

最好， 强尼

Answer 2

这不能回答您的问题，但它会给出一个您可以在网页中进行交互的结果：不要使用plot_ly，使用rgl。例如，

library(rgl)

# Your initial values...

r <- function(x) x^2
int <- c(1, 3)
nx <- 20
ntheta <- 36

# Set up x and colours for each x

x <- seq(int[1], int[2], length.out = nx)
cols <- colorRampPalette(c("blue", "yellow"), space = "Lab")(nx)

clear3d()
shade3d(turn3d(x, r(x), n = ntheta,  smooth = TRUE, 
        material = list(color = rep(cols, each = 4*ntheta))))
aspect3d(1,1,1)
decorate3d()
rglwidget()

您可以通过一些小技巧在颜色上做得更好：您可能想要创建一个使用x 或r(x) 来设置颜色的函数，而不是像我一样只是重复颜色。

结果如下：

Answer 3

我有另一个破解，并有一个更接近的解决方案，使用“表面”类型。有帮助的是查看您的第一个曲面图的结果，nx = 5 和 ntheta = 18。它的原因是因为它连接 zs 中的列（跨 x 点）的方式。它必须从它周围的较大环的一部分向上连接，这会导致密度飙升以满足这一点。

我无法 100% 摆脱这种粗鲁的行为。我进行了以下更改：

1. 在边缘周围的 theta 上添加一些小点：两个密度连接的地方。这减小了锯齿部分的大小，因为边界附近还有更多点
2. 从 zs 到 zs2 的计算：通过添加 0 来确保每个环与外部环的维度相同。
3. 将 nx 增加到 40 并将 ntheta 减少到 18 - 更多的 x 使步长更小。减少运行时间的 ntheta，因为我已经添加了更多点

步骤在于它如何尝试连接 x 环。理论上，如果你有更多的 x 环，它应该可以消除这种粗糙，但运行起来很耗时。

我认为这不能 100% 回答问题，我不确定这个库是否最适合这项工作。如果有任何问题，请与我们联系。

library(dplyr)
library(plotly)

# radius depends on x
r <- function(x) x^2

# interval of interest
int <- c(1, 3)

# number of points along the x-axis
nx <- 40

# number of points along the rotation
ntheta <- 18

# set x points and get corresponding radii
coords <- data_frame(x = seq(int[1], int[2], length.out = nx), r = r(x))

# theta: add small increments at the extremities for the density plot
theta <- seq(0, pi, length.out = ntheta / 2 + 1)
theta <- c(theta, pi + theta)
theta <- theta[theta != 2*pi]
inc <- 0.00001
theta <- c(theta, inc, pi + inc, pi - inc, 2*pi - inc)
theta <- sort(theta)

coords %<>%
  rowwise() %>%
  do(data_frame(x = .$x, r = .$r, theta = theta)) %>%
  ungroup %>%
  mutate(y = r * cos(theta), z = r * sin(theta))

# get all x & y values used (sort to connect halves on the side)
xs <-
  unique(coords$x) %>%
  sort
ys <-
  unique(coords$y) %>%
  sort

# for each possible x/y pair: get z^2 value
coords <-
  expand.grid(x = xs, y = ys) %>%
  as_data_frame %>%
  mutate(r = r(x), z2 = r^2 - y^2)

# format z coordinates above x/y plane as matrix where columns
# represent x and rows y
zs <- matrix(sqrt(coords$z2), ncol = length(xs), byrow = TRUE)
zs2 <- zs

L <- ncol(zs)
for(i in (L-1):1){
  w <- which(!is.na(zs[, (i+1)]) & is.na(zs[, i]))
  zs2[w, i] <- 0
}

# format x coordiantes as matrix as above (for color gradient)
gradient <-
  rep(xs, length(ys)) %>%
  matrix(ncol = length(xs), byrow = TRUE)

# plot upper half of shape as surface
p <- plot_ly(x = xs, y = ys, z = zs2, surfacecolor = gradient,
             type = "surface", colorbar = list(title = 'x'))

# plot lower have of shape as second surface
p %>%
  add_surface(z = -zs2, showscale = FALSE)

Answer 4

一种解决方案是翻转您的轴，以便您围绕 z 轴而不是 x 轴旋转。我不知道这是否可行，考虑到您要添加这个数字的现有图表，但它确实很容易解决“牙齿”问题。

xs <- seq(-9,9,length.out = 20)
ys <- seq(-9,9,length.out = 20)

coords <-
  expand.grid(x = xs, y = ys) %>%
  mutate(z2 = (x^2 + y^2)^(1/4))

zs <- matrix(coords$z2, ncol = length(xs), byrow = TRUE)

plot_ly(x = xs, y = ys, z = zs, surfacecolor = zs,
             type = "surface", colorbar = list(title = 'x')) %>% 
  layout(scene = list(zaxis = list(range = c(1,3))))