在 R 或 Matlab 中绘制这样颜色的多边形

Question

http://www.texample.net/tikz/examples/lindenmayer-systems/

我的示例代码如下所示，我不知道如何用色相着色。

plot.koch <- function(k,col="blue"){ 
  plot(0,0,xlim=c(0,1), ylim=c(-sqrt(3)/6,sqrt(3)/2), asp = 1,type="n",xlab="", ylab="")
  plotkoch <- function(x1,y1,x2,y2,n){
    if (n > 1){
      plotkoch(x1,y1,(2*x1+x2)/3,(2*y1+y2)/3,n-1); 
      plotkoch((2*x1+x2)/3,(2*y1+y2)/3,(x1+x2)/2-(y1-y2)*sqrt(3)/6,(y1+y2)/2-(x2-x1) *sqrt(3)/6,n-1);
      plotkoch((x1+x2)/2-(y1-y2)*sqrt(3)/6,(y1+y2)/2-(x2-x1)*sqrt(3)/6,(2*x2+x1)/3,(2 *y2+y1)/3,n-1); 
      plotkoch((2*x2+x1)/3,(2*y2+y1)/3,x2,y2,n-1)
    }    
    else { 
      x=c(x1,(2*x1+x2)/3,(x1+x2)/2-(y1-y2)*sqrt(3)/6,(2*x2+x1)/3,x2); 
      y=c(y1,(2*y1+y2)/3,(y1+y2)/2-(x2-x1)*sqrt(3)/6,(2*y2+y1)/3,y2); 
      polygon(x,y,type="l",col=col) 
    }
  }
  plotkoch(0,0,1,0,k) 
  plotkoch(0.5,sqrt(3)/2,0,0,k) 
  plotkoch(1,0,0.5,sqrt(3)/2,k)
}
  plot.koch(3, col=3)

Answer 1

我会这样做：

对于任何绘制的像素获取其位置x,y
计算 angle=atan2(y-y0,x-x0)

其中 x0,y0 是科赫雪花的中间位置
根据角度计算颜色

如果你使用 HSV 然后 hue=angle 并计算目标颜色值（我假设 RGB）。如果你想要可见光谱颜色，你可以试试我的：
- RGB values of visible spectrum
只需将角度范围angle=<0,2*Pi> [rad]转换为波长l=<400,700> [nm]即可：
```
l = 400.0 + (700.0-400.0)*angle/(2.0*M_PI)
```
渲染像素

[备注]

没有使用 R 也没有使用 Matlab 所以你需要自己编写代码。角度可能需要一些移动以匹配您的坐标系，例如：

const angle0=???; // some shift constant [rad]
angle+=angle0; // or angle=angle0-angle; if the direction is oposite
if (angle>=2.0*M_PI) angle-=2.0*M_PI;
if (angle<      0.0) angle+=2.0*M_PI;

如果您将其绘制为多边形，那么您需要计算每个顶点的颜色而不是每个像素，但是您可能会遇到问题，因为这不是凸多边形。那么如何保证中点颜色？？？恐怕你需要使用某种三角测量，因为简单的三角扇会失败......

唯一明显的是填充整个颜色space，然后用黑色绘制轮廓，然后用白色从外部填充所有非黑色像素

Answer 2

这是我尝试解决您的问题的方法。目前它也在雪花之外绘制颜色。如果你能弄清楚点是在雪花内部还是外部，你应该能够只删除 df_fill 中的外部点。在这里，我首先创建用于绘制多边形的 data.frame。然后我为背景颜色创建 data.frame。最后我使用 ggplot2 来绘制数据。

# creating relevant data
data.koch <- function(k){ 
  df <- data.frame(x = 0, 
                   y = 0, 
                   grp = 0)
  plotkoch <- function(x1, y1, x2, y2, n, data){
    if (n==1) {
      x=c(x1,(2*x1+x2)/3,(x1+x2)/2-(y1-y2)*sqrt(3)/6,(2*x2+x1)/3,x2) 
      y=c(y1,(2*y1+y2)/3,(y1+y2)/2-(x2-x1)*sqrt(3)/6,(2*y2+y1)/3,y2) 
      df <- rbind(data, data.frame(x, y, grp=max(data$grp)+1))
    }
    if (n > 1){
      df <- plotkoch(x1,y1,(2*x1+x2)/3,(2*y1+y2)/3,n-1, data = data)
      df <- plotkoch((2*x1+x2)/3,(2*y1+y2)/3,(x1+x2)/2-(y1-y2)*sqrt(3)/6,(y1+y2)/2-(x2-x1) *sqrt(3)/6,n-1, data=df)
      df <- plotkoch((x1+x2)/2-(y1-y2)*sqrt(3)/6,(y1+y2)/2-(x2-x1)*sqrt(3)/6,(2*x2+x1)/3,(2 *y2+y1)/3,n-1, data=df) 
      df <- plotkoch((2*x2+x1)/3,(2*y2+y1)/3,x2,y2,n-1, data=df)
    }    
    return(df)
  }
  df <- plotkoch(0,0,1,0,k, data = df) 
  df <- plotkoch(0.5,sqrt(3)/2,0,0,k, data = df) 
  df <- plotkoch(1,0,0.5,sqrt(3)/2,k, data = df)
  return(df)
}
# plotting functon
plot.koch <- function(k){
  stopifnot(require(ggplot2))
  if (is.data.frame(k)) df <- k
  else df <- data.koch(k)
  # filling data (CHANGE HERE TO GET ONLY INSIDE POINTS)
  l <- 500
  df_fill <- expand.grid(x=seq(0, 1, length=l), 
                         y=seq(-sqrt(3)/6, sqrt(3)/2, length=l))
  df_fill[, "z"] <- atan2(-df_fill[, "y"] + sqrt(3)/6, df_fill[, "x"] - 0.5) + pi/2
  df_fill[df_fill[, "z"] < 0, "z"] <- df_fill[df_fill[, "z"] < 0, "z"] + 2*pi
  # plotting
  ggplot(df, aes(x, y, group=grp)) + 
    geom_raster(data = df_fill, 
                aes(fill=z, group=NULL), 
                hjust = 0,
                vjust = 0,
                linetype='blank') + 
    geom_path(data=df, size=1) +
    scale_fill_gradientn(colours = rainbow(30), guide = 'none') + 
    scale_x_continuous(name = '', limits = c(0, 1), expand=c(0, 0)) + 
    scale_y_continuous(name = '', limits = c(-sqrt(3)/6,sqrt(3)/2), expand=c(0, 0)) +
    coord_fixed() +
    theme_bw() +
    theme(axis.line = element_blank(), 
          panel.grid = element_blank(), 
          axis.ticks = element_blank(), 
          axis.text = element_blank())
}
#
p <- plot.koch(4)
print(p)

Answer 3

这是一种在 R 中使用空间对象的方法，混合了 sp、rgeos 和 raster 包。

对 return 函数的轻微修改 x，y 坐标给用户（并以正确的顺序）：

koch <- function(k) { 
  yy <- xx <- numeric(0)
  Koch <- function(x1, y1, x2, y2, n) {
    if (n > 1){
      Koch(x1, y1, (2*x1+x2)/3, (2*y1+y2)/3, n-1); 
      Koch((2*x1+x2)/3, (2*y1+y2)/3, (x1+x2)/2-(y1-y2)*sqrt(3)/6, 
           (y1+y2)/2-(x2-x1) *sqrt(3)/6, n-1);
      Koch((x1+x2)/2-(y1-y2)*sqrt(3)/6, (y1+y2)/2-(x2-x1)*sqrt(3)/6, 
           (2*x2+x1)/3, (2 *y2+y1)/3, n-1); 
      Koch((2*x2+x1)/3, (2*y2+y1)/3, x2, y2, n-1)
    }    
    else { 
      x <- c(x1, (2*x1+x2)/3, (x1+x2)/2-(y1-y2)*sqrt(3)/6, (2*x2+x1)/3, x2); 
      xx <<- c(xx, x)
      y <- c(y1, (2*y1+y2)/3, (y1+y2)/2-(x2-x1)*sqrt(3)/6, (2*y2+y1)/3, y2); 
      yy <<- c(yy, y)
     }
  }
  Koch(0, 0, 1, 0, k)
  Koch(1, 0, 0.5, sqrt(3)/2, k)
  Koch(0.5, sqrt(3)/2, 0, 0, k) 
  xy <- data.frame(x=xx, y=yy)
  rbind(unique(xy), xy[1, ])
}

创建色带：

colr <- colorRampPalette(hcl(h=seq(0, 360, len=100), c=100))

使用koch函数获取顶点：
```
xy <- koch(4)
```

加载空间包并从分形创建SpatialPolygons对象并绘制一次以设置绘图区域。

library(sp)
library(rgeos)
library(raster)

poly <- SpatialPolygons(list(Polygons(list(Polygon(xy)), 1)))
plot(poly)

绘制一系列具有所需原点和足够大半径以覆盖分形多边形的线段（此处我们使用半径 r <- 1）。

r <- 1
mapply(function(theta, col) {
  segments(0.5, 0.3, 0.5 + r*cos(theta), 0.3 + r*sin(theta), lwd=3, col=col) 
}, seq(0, 360, length=1000)*pi/180, colr(1000))

创建绘图区域和分形多边形之间差异的第二个多边形，并绘制它（使用 col='white'）以掩盖不需要的渐变区域。
```
plot(gDifference(as(extent(par('usr')), 'SpatialPolygons'), poly), 
     col='white', border='white', add=TRUE)
```
再次绘制多边形。
```
plot(poly, add=TRUE)
```

Answer 4

这是我用网格包的解决方案。

##data
  koch <- function(k) {
        yy <- xx <- numeric(0)
        Koch <- function(x1, y1, x2, y2, n) {
            if (n > 1) {
                Koch(x1, y1, (2 * x1 + x2)/3, (2 * y1 + y2)/3, n - 1)
                Koch((2 * x1 + x2)/3, (2 * y1 + y2)/3, (x1 + x2)/2 - (y1 - 
                  y2) * sqrt(3)/6, (y1 + y2)/2 - (x2 - x1) * sqrt(3)/6, 
                  n - 1)
                Koch((x1 + x2)/2 - (y1 - y2) * sqrt(3)/6, (y1 + y2)/2 - 
                  (x2 - x1) * sqrt(3)/6, (2 * x2 + x1)/3, (2 * y2 + y1)/3, 
                  n - 1)
                Koch((2 * x2 + x1)/3, (2 * y2 + y1)/3, x2, y2, n - 1)
            } else {
                x <- c(x1, (2 * x1 + x2)/3, (x1 + x2)/2 - (y1 - y2) * sqrt(3)/6, 
                  (2 * x2 + x1)/3, x2)
                xx <<- c(xx, x)
                y <- c(y1, (2 * y1 + y2)/3, (y1 + y2)/2 - (x2 - x1) * sqrt(3)/6, 
                  (2 * y2 + y1)/3, y2)
                yy <<- c(yy, y)
            }
        }
        Koch(0, 0, 1, 0, k)
        Koch(1, 0, 0.5, sqrt(3)/2, k)
        Koch(0.5, sqrt(3)/2, 0, 0, k)
        xy <- data.frame(x = (xx - min(xx))/(max(xx) - min(xx)), y = (yy - 
            min(yy))/(max(yy) - min(yy)))
        rbind(unique(xy), xy[1, ])
    }
xy <- koch(5)
##Plot
library(grid)
grid.newpage()
pushViewport(dataViewport(xy$x, xy$y), plotViewport(c(1, 1, 1, 1)))
    for (i in 1:nrow(xy)) {
       grid.path(x = c(xy[i, 1], xy[i + 1, 1], mean(xy$x)), 
                 y = c(xy[i, 2], xy[i + 1, 2], mean(xy$y)), 
                 gp = gpar(col = rainbow(nrow(xy))[i], 
                           fill = rainbow(nrow(xy))[i]))
        }

在 R 或 Matlab 中绘制这样颜色的多边形

Draw a polygon colored like this in R or Matlab

graphics

matlab

plot

r

polygons