在R中绘制光栅图像，画一个圆并遮盖圆外的像素

Question

下面的代码绘制了一张图片，然后在该图片上绘制了一个圆。我想把落在圆外的所有像素都变成黑色。我该如何实现？

library(raster)
library(plotrix)
r1 <- brick(system.file("external/rlogo.grd", package="raster"))
width=50
height=40
x <- crop(r1, extent(0,width,0,height))
plotRGB(x)
circlex=20
circley=15
radius=10
draw.circle(circlex,circley,radius,border="blue")

- user2543622

Answer 1

使用str()查看'x'对象，你会看到这个：

..@ data    :Formal class '.MultipleRasterData' [package "raster"] with 14 slots
  .. .. ..@ values    : num [1:2500, 1:3] 255 248 221 199 198 210 221 190 104 79 ...
  .. .. .. ..- attr(*, "dimnames")=List of 2
  .. .. .. .. ..$ : NULL
  .. .. .. .. ..$ : chr [1:3] "red" "green" "blue"

...所以1:50乘1:50的值被映射到三列。X值可能是0:2500 %% 50，y值可能是0:2500 %/% 50

因此请记住，光栅对象的“原点”为左上角，而绘图函数的“原点”为左下角，因此y值为20变为50-20或30，这就给出了你所要求的结果（抱歉将y序列放在第一位）：

x@data@values[( ((1:2500 %/% 50 )- 30)^2 + ((1:2500 %% 50) - 20)^2 ) >=100, 1] <- 0
 x@data@values[( ((1:2500 %/% 50 )- 30)^2 + ((1:2500 %% 50) - 20)^2 ) >=100, 2] <- 0
 x@data@values[( ((1:2500 %/% 50 )- 30)^2 + ((1:2500 %% 50) - 20)^2 ) >=100, 3] <- 0
 plotRGB(x)
 draw.circle(20,20,10,border="blue")

The logic is that the criteria are of the form (x-dx)^2+(y-dy)^2 > r^2 where dx and dy are the center coordinates of the circle and r is the radius == 10.

对于每个颜色层，具有命名参数的代码将类似于制作最暗的“红色”代码。这会给出一个大致的圆形掩模，尽管正确处理中心线不是很好。

x@data@values[( ((1:(height*width) %/% (height*5/4) )- (height-circley*5/4) )^2 + 
            ((1:(height*width) %% width) -       circlex )^2 ) >= radius^2, 1] <- 0

进一步的实验得出了这个结果，看起来非常接近：

x@data@values[( ((1:(height*width) %/% (height) )- (height-circley) *5/4)^2 + 
                 ((1:(height*width) %% width) -       circlex )^2 ) >= radius^2, 1] <- 0
plotRGB(x, asp=5/4, addfun=function() draw.circle(circlex,circley,radius,border="blue") )

显然，你可以在所有出现5/4的地方将宽度/高度缩放因子替换为新的纵横比。

- IRTFM

---

你是那个没有命名参数的问题提出者。 - IRTFM

可能需要引入一个缩放因子，因为“像素”的宽高比不再是1.0。 - IRTFM

我再次检查了一下。@file@nbands的值是3，所以我认为图像仍然有3个波段。@data@values的值是logi[0,0]。因此，我很难获取颜色值。我正在使用r1 <- brick("/mypath/results_circle.jpg")读取JPEG图像。 - user2543622

你不知怎么地在@data@values中得到了一个空的逻辑矩阵。如何发生这种情况是一个需要解决的难题。你可能需要在另一个问题中发布str(r1)的完整输出。我真的认为我已经完全回答了这个问题。 - IRTFM

如果我知道那是怎么发生的，我就会回答了。 - IRTFM

Answer 2

这是一个更通用的解决方案，也适用于具有不同坐标系统的栅格。

函数 rasterToPoints() 将栅格转换为点。按照您的示例，head(rasterToPoints(x)) 返回以下内容：

> head(rasterToPoints(x))
       x    y red green blue
[1,] 0.5 39.5 255   255  251
[2,] 1.5 39.5 204   205  199
[3,] 2.5 39.5 171   172  164
[4,] 3.5 39.5 157   159  148
[5,] 4.5 39.5 162   164  151
[6,] 5.5 39.5 187   189  176

我们需要找出哪些点落在圆外，并将它们的值设为零：

is_outside_circ = (rasterToPoints(x)[,1] - circlex)^2 + (rasterToPoints(x)[,2] - circley)^2 >= radius^2
x@data@values[ is_outside_circ,] <- 0
plotRGB(x)
draw.circle(circlex,circley,radius,border="blue")

结果：圆外黑点

- jsga