0

我有多年的栅格列表（.tif 格式）。这是来自 landsat 的 16 天 NDVI，我想制作每月 NDVI（两个连续栅格的平均值）并将其作为每月平均值保存在相同或不同的目录中

我已经列出了栅格栅格并将其堆叠起来，后来我使用 stackApply 来计算平均值，但它会产生空栅格。我一年有 23 张图像，我想平均它并制作 12 个月。这就是我的光栅文件的样子

 "landsatNDVISC05SLC2000001.tif" "landsatNDVISC05SLC2000017.tif"
 "landsatNDVISC05SLC2000033.tif" "landsatNDVISC05SLC2000049.tif"
 "landsatNDVISC05SLC2000065.tif" "landsatNDVISC05SLC2000081.tif"
 "landsatNDVISC05SLC2000097.tif" "landsatNDVISC05SLC2000113.tif"
 "landsatNDVISC05SLC2000129.tif" "landsatNDVISC05SLC2000145.tif"
 "landsatNDVISC05SLC2000161.tif" "landsatNDVISC05SLC2000177.tif"
 "landsatNDVISC05SLC2000193.tif" "landsatNDVISC05SLC2000209.tif"
 "landsatNDVISC05SLC2000225.tif" "landsatNDVISC05SLC2000241.tif"
 "landsatNDVISC05SLC2000257.tif" "landsatNDVISC05SLC2000273.tif"
 "landsatNDVISC05SLC2000289.tif" "landsatNDVISC05SLC2000305.tif"
 "landsatNDVISC05SLC2000321.tif" "landsatNDVISC05SLC2000337.tif"
 "landsatNDVISC05SLC2000353.tif

此代码有效，但会产生超过十二个空栅格，我还想将栅格砖保存为单个子集每月栅格

library(raster)
lrast<-list.files("G:/LANDSAT-NDVI/testAverage")
layers<-paste("landsatNDVISC05SLC2000", seq(from=001, to=353,by=16))
stak<-stack(lrast)
raster<-stackApply(stak, layers, fun = mean)

我想将 landsatNDVISC05SLC2000001.tif 和 landsatNDVISC05SLC2000017.tif 作为 landsatNDVISC05SLC2000M1.tif 的月平均值。同样，33,49 并且由于我只有 23 个栅格，我想将 landsatNDVISC05SLC2000353.tif 保留为 landsatNDVISC05SLC2000M12.tif

块引用

4

1 回答 1

0

不确定 stackapply 是如何工作的，但是这样的东西应该可以做需要的东西。

library(raster)
files <- list.files(path = "...", full.names = T, pattern = ".tif")

stk <- stack()

for (i in files){
  print(i)
  as <- raster(files[i])
  stk <- addLayer(stk, as)
}

jday <-c("landsatNDVISC05SLC2000017.tif","landsatNDVISC05SLC2000033.tif",
"landsatNDVISC05SLC2000049.tif","landsatNDVISC05SLC2000065.tif","landsatNDVISC05SLC2000081.tif",
"landsatNDVISC05SLC2000097.tif","landsatNDVISC05SLC2000113.tif","landsatNDVISC05SLC2000129.tif",
"landsatNDVISC05SLC2000145.tif","landsatNDVISC05SLC2000161.tif","landsatNDVISC05SLC2000177.tif",
"landsatNDVISC05SLC2000193.tif","landsatNDVISC05SLC2000209.tif","landsatNDVISC05SLC2000225.tif",
"landsatNDVISC05SLC2000241.tif","landsatNDVISC05SLC2000257.tif","landsatNDVISC05SLC2000273.tif",
"landsatNDVISC05SLC2000289.tif","landsatNDVISC05SLC2000305.tif","landsatNDVISC05SLC2000321.tif",
"landsatNDVISC05SLC2000337.tif","landsatNDVISC05SLC2000353.tif")

jday <- as.numeric(substr(jday, 24, 25)) #substract the julien days (which I think these number represent before .tif; or you can substract the names from the 'files' vector)

dates <- as.Date(jday, origin=as.Date("2000-01-01")) # create a Date vector

stk <- setZ(stk, dates) # assign the date vector to the raster stack

raster <- zApply(stk, by = format(dates,"%Y-%m"), fun = mean, na.rm = T) # create the monthly stack
于 2019-07-05T08:34:30.093 回答