Как нарисовать трехмерный граф в R

Question

f(x,y)= (1/25)*(20-x)/x   10<x<20, x/2 <y <x
        0                   o.t

Я должен создать это изображение с помощью этого выражения.

но

x <- seq(10, 20, length=20)
y <- seq(10, 20, length=20)
f <- function(x,y){(1/25)*(20-x)/5}
z <- outer(x,y,f)
persp(x,y,z,theta=30,phi=30, expand=0.5,col=rainbow(19), border=NA)

что не так?

Answer 1

@ Ответ SRhm, вероятно, является лучшим выбором, но если вы хотите жить на переднем крае, вы можете избавиться от зазубренного диагонального края, используя версию разработки rgl (из R-forge), по крайней мере, версию0.100.8.

Эта версия поддерживает триангуляции с границами с использованием пакета tripack.Таким образом, вы устанавливаете сетку значений в диапазоне xy, затем определяете границы области, используя уравнения, и вы получаете гладкие края.Например:

library(tripack)
library(rgl)
g <- expand.grid(x=10:20, y=5:20)
keep <- with(g, 10 < x & x < 20 & x/2 < y & y < x)
g2 <- g[keep,]
tri <- tri.mesh(g2)

# Set up boundary constraints
cx <- c(10:20, 20: 10)
cy <- c(seq(5, 10, len=11), 20:10)
tri2 <- add.constraint(tri, cx, cy, reverse = TRUE)

# This isn't necessary, but shows where the formula will be evaluated
plot(tri2)

Возможно, было бы лучше заполнить некоторые из левого и правого краев большим количеством точек, чтобы избежать этих больших треугольников,но пропустите это сейчас.

z <- with(tri2, (1/25)*(20-x)/x)
# Now plot it, using the map2color function @SRhm found:
#source: /10119043/plot3d-kak-izmenit-tsvet-poverhnosti-osi-z-na-tsvet-teplovoi-karty
map2color <- function(x, pal, limits = range(x,na.rm=T)){
  pal[findInterval(x, seq(limits[1], limits[2], length.out = length(pal) + 1), 
               all.inside=TRUE)]
}
persp3d(tri2, z, col = map2color(z, rainbow(100)))

После поворота вы получите это представление:

Answer 2

Вы должны замаскировать z на основе ограничения.В качестве предложения вы можете использовать удивительный интерактивный пакет rgl в R.

#source: /10119043/plot3d-kak-izmenit-tsvet-poverhnosti-osi-z-na-tsvet-teplovoi-karty
map2color <- function(x, pal, limits = range(x,na.rm=T)){
  pal[findInterval(x, seq(limits[1], limits[2], length.out = length(pal) + 1), 
               all.inside=TRUE)]
}

x <- seq(10, 20, length=20)
y <- seq(10, 20, length=20)
mask <- sapply(x,function(m) sapply(y,function(n) if((n>m/2)&(n<m)){(1/25)*(20-m)/5}else{ NA }))
z <- outer(x,y,f)
z <- z * mask
#persp(x,y,z, col= map2color(z, rainbow(100)),border = NA)
library(rgl)
persp3d(x,y,z,col = map2color(z, rainbow(100)),theta=30,phi=30)