R对带有正方形多边形的栅格进行采样

Question

我想通过创建小栅格 100x100 单元格来采样大栅格。 我不知道该怎么做，所以欢迎任何想法

我的实际线索：

library(raster)
library(spatstat)
library(polyCub)

r <- raster(ncol=1000,nrow=1000) # create empty raster
r[] <- 1:(1000*1000)             # Raster for testing
e <- extent(r)                   # get extend
# coerce to a SpatialPolygons object
p <- as(e, 'SpatialPolygons')  


nc <- as.owin.SpatialPolygons(p) #polyCub
pts <- rpoint(50, win = nc)
plot(pts)

现在我需要在我的 50 个点周围生成 100x100 单元格正方形，我想使用这些正方形裁剪 r 并单独堆叠每个小栅格...

Answer 1

@adrian-baddeley 的回答基本上包含了做什么 你要。如果您只是想要一个包含小 im 对象的列表 您只需将 im 对象通过 owin 对象子集到 100x100 框 提取相关区域。这是一个例子（点数较少 避免过度绘制）

library(raster)
library(spatstat)
library(maptools)

r <- raster(ncol=1000,nrow=1000) # create empty raster
r[] <- 1:(1000*1000)             # Raster for testing
e <- extent(r)                   # get extend
# coerce to a SpatialPolygons object
p <- as(e, 'SpatialPolygons')  

nc <- as.owin.SpatialPolygons(p)
set.seed(42)
pts <- rpoint(7, win = nc)

rim <- as.im.RasterLayer(r)
Box <- owin(c(-50,50) * rim$xstep, c(-50,50) * rim$ystep)

以下是大小为 100x100 的 im 对象列表

imlist <- solapply(seq_len(npoints(pts)),
                   function(i) rim[shift(Box, pts[i])])

这是该区域中im 对象和顶部点的图

plot(pts)
for(i in imlist) plot(i, add = TRUE)
plot(pts, pch = 19, add = TRUE)

您可以使用

转换为栅格图层列表

rasterList <- lapply(imlist, as, Class = "RasterLayer")

PS：以下是im对象的原始大小列表 NA 如果您需要该格式，请在 100x100 框外

imlist <- solapply(seq_len(npoints(pts)),
                   function(i) rim[shift(Box, pts[i]), drop = FALSE])

Answer 2

如果要使用spatstat，则需要将光栅对象r 转换为spatstat 支持的im 类对象。您可以在 maptools 包中进行此转换。将此图像对象称为rim。然后你可以这样做

Box <- owin(c(-50,50) * rim$xstep, c(-50,50) * rim$ystep)
BoxesUnion <- MinkowskiSum(pts, Box)
W <- intersect.owin(as.mask(rim), BoxesUnion)

这将为您提供被正方形覆盖的栅格子集。 如果您想保持正方形分开，请执行以下操作

M <- as.mask(rim)
BoxList <- solapply(seq_len(npoints(pts)), 
                      function(i) intersect.owin(M, shift(Box, pts[i])))

那么BoxList 是各个子栅格的列表。