如何使用 R 将 shapefile 剪辑到我的光栅？

Question

我试图将光栅屏蔽到形状文件边界，但出现错误。我怎样才能正确地执行这个面具？

原始数据可以在这里找到，标题为“data_for_question.txt”。它经过格式化，以便用户可以（从 Web 应用程序）将文本直接复制并粘贴到 R 窗口中并生成数据框。否则，如果不想生成数据，也可以在同一链接中找到输出栅格 (example_raster.tif) 和 shapefile (field_boundary.shp)。

这是我尝试过的：

#Import necessary libraries
library(pacman)
p_load(sf,
       spatstat,
       maptools,
       tidyverse,
       ggplot2, 
       gstat,
       sp,
       rgdal,
       raster,
       spdep)

#Read shapefile
shp <- st_read("field_boundary.shp")

#Generate data to run interpolation on and project it to the desired CRS
data_sp <- SpatialPointsDataFrame(coords, 
                                      data[, c("OM", "data2")], 
                                      proj4string = CRS('+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs'))


#Perform an IDW interpolation:
grd <- SpatialPixels(SpatialPoints(makegrid(data_sp, n=10000)), proj4string = proj4string(data_sp)) #Generate grid for interpolation

plot(grd)

interp <- idw(formula = OM ~ 1, data_sp, grd, idp = 0.5, nmax = 12)
plot(interp) #Makes for a very pretty picture!

#Convert to raster
rast <- raster(interp)
plot(rast)
shp <- st_transform(shp, crs(rast))

#Crop and mask the raster
crop_rast <- crop(rast, shp)
crop_om <- mask(crop_rast, mask = shp)

错误发生在这里：

Error in h(simpleError(msg, call)) : 
  error in evaluating the argument 'x' in selecting a method for function 'addAttrToGeom': sp supports Z dimension only for POINT and MULTIPOINT.
use `st_zm(...)` to coerce to XY dimensions

Answer 1

错误消息表明我们需要将几何图形强制为 XY 维度。因此，在使用 shapefile 裁剪和屏蔽栅格之前，请尝试执行此操作：

shp<- st_zm(shp)

然后继续裁剪和遮罩光栅。