使用多边形裁剪栅格文件并使用相同的文件名写入输出

Question

我有很多光栅文件（卫星图像，都可以在 geotiff .tif 扩展名中找到）。有些文件作为单个文件被分成所有波段，有些文件有多个波段。由于这在我的硬盘上使用了很多space，我想用我感兴趣的区域裁剪每个文件，我有一个 shapefile 多边形。

我接近我自己的解决方案并使用以下代码将裁剪后的图像作为新的 .tif 文件获取：

library(raster)
rasterfiles = list.files(path=getwd(), pattern = "*.TIF", full.names=TRUE)
s = stack(rasterfiles)
shp = readOGR("Area.shp")
rasterfiles_crop = crop(s, extent(shp))

output = writeRaster(rc, 'out.tif', format="GTiff", overwrite=TRUE, bylayer = TRUE)

使用此代码，我收到文件名 out_1.tif、out_2.tif 等...

不幸的是，生成的文件只有 1 个波段，因此 R 仅在涉及多波段 TIF 图像时识别第一个波段。

我想保留所有波段和原始文件名，只是在新波段的末尾添加“_crop”。有人可以在这里帮助我如何更改代码吗？

谢谢

Answer 1

你可以把它们写成一个循环

library(raster)
library(rgdal)

shp <- readOGR("Area.shp")     

infiles <- list.files(path=getwd(), 
                      pattern="*.TIF", 
                      full.names=TRUE)

outfiles <- file.path(YourOutputPath, 
                      paste0(basename(tools::file_path_sans_ext(infiles)),
                      "_crop.tif")
                      )

for (i in seq_along(infiles)) {
  r <- crop(raster(infiles[i]), shp)
  writeRaster(r, filename=outfiles[i])
}

Answer 2

我现在找到了一个解决方案，下面的代码列出了一个文件夹中的所有 TIF 文件，并且多波段 tif 在裁剪过程后保留其波段：

library(raster)
library(rgdal)

setwd("input-folder")

 ## polygon with crop-extend ##
shp <- readOGR("area.shp")

 ## load tif files ##
infiles = list.files(path=getwd(),
                     pattern="*.tif$|*.TIF$")

 ## Filenames with desired suffix and output place ##
outfiles = file.path("D:/Downloads/BDA/Output",
                     paste0(basename(tools::file_path_sans_ext(infiles)),
                     ".tif"))

 ## crop and output settings (compression and datatype)
 for (i in seq_along(infiles)) {
      r = crop(stack(infiles[i]), shp)
      writeRaster(r, filename=outfiles[i],
                  bylayer=FALSE,
                  format="GTiff",
                  datatype="INT1U",
                  options="COMPRESS=ZIP",
                  overwrite=TRUE)
 }

感谢 Richard 提供的循环代码！

关于数据类型：如果 R 可以检查输入文件的数据类型并自动为裁剪输出选择相同的数据类型，那就太好了。现在我必须手动指定数据类型。否则即使输入文件只有 8 位无符号 (INT1U) 或 16 位有符号 (INT2S)，输出文件也将是 float32 (FLT4S)。

 datatype= same.as.input.file

使用多边形裁剪栅格文件并使用相同的文件名写入输出

Crop raster files with Polygon and writing the output with the same filename

crop

polygon

raster

geotiff