У меня есть данные MODIS, которые содержат продукт за 8 дней. Некоторые месяцы имеют 3 изображения, а некоторые месяцы имеют четыре изображения. Итак, как я могу создать сценарий R для расчета ежемесячной суммы из этих данных?
Я использую набор данных MODIS MOD17A2HGF GPP, который представляет собой сетку 500x500 метров с суммой GPP за 8 дней. Первый файл заканчивается на _001, что означает, что он начинается 1 января и продолжается до 8 числа. Следующий файл, _009, охватывает период с 8 по 16 января и так далее. введите описание изображения здесь
Теперь я хочу суммировать GPP по месяцам (с января по декабрь).
В этой папке 920 файлов (20 лет), из которых каждый год 46 файлов.
В ArcGIS я знаю, как суммировать с помощью растрового калькулятора. Но если я буду суммировать каждый месяц, он будет делать столько повторяющейся работы.
введите описание изображения здесь
В R я пытался сделать это именно так
library(terra)
a <- rast('MOD17A2HGF_GPP_2001_001.tif')
b <- rast('MOD17A2HGF_GPP_2001_009.tif')
c <- rast('MOD17A2HGF_GPP_2001_017.tif')
d <- rast('MOD17A2HGF_GPP_2001_025.tif')
GPP_January <- a+b+c+d
writeRaster(GPP_January,'F:/test/GPP_January_2001.tif')
Он также будет выполнять так много повторяющихся работ, если я буду суммировать каждый месяц. Я не знаю, как выбрать указанные файлы месяца в цикле R.
Я знаю, как суммировать годовой GPP
library(tidyverse)
library(terra)
GPP2001 <- list.files(pattern = '2001',full.names = T) %>%
rast() %>% sum()
Но как я могу суммировать ежемесячные GPP в год в R или ArcGIS?
Получить все файлы/слои
library(terra)
ff <- list.files(pattern = "\\.tif$", full.names = T)
r <- rast(ff)
Для каждого файла извлеките год и номер дня.
# ff <- c('MOD17A2HGF_GPP_2001_001.tif', 'MOD17A2HGF_GPP_2001_009.tif')
s <- gsub('MOD17A2HGF_GPP_', "", basename(ff))
s <- gsub('.tif', "", s)
s <- strsplit(s, "_")
s <- do.call(rbind, s)
s
# [,1] [,2]
#[1,] "2001" "001"
#[2,] "2001" "009"
Получить месяц и объединить с годом
d <- as.Date(as.numeric(s[,2]), origin = paste(s[,1], "-01-01", sep = "")) - 1
m <- format(d, "%m")
ym <- paste0(s[,1], m)
ym
#[1] "200101" "200101"
Теперь вы можете рассчитать среднемесячное значение следующим образом.
x <- tapp(x, ym, mean)
(точнее, среднее значение для дат каждого месяца; это не совсем среднемесячное значение, но, вероятно, достаточно хорошее).
Я не могу описать, насколько я взволнован в данный момент. Меня вдохновил аналогичный вопрос, на который вы ответили в 2018 году, и я подумал, что это идеально, но этот вопрос касался средних значений, и я не знаю, как изменить последнюю строку кода. Я не ожидал ответа оригинального автора вскоре после моего вопроса, который так волнует!
Я использую mean
(не sum
), потому что вы заявили, что хотите вычислить «среднемесячное значение».
Я сожалею об этом, это ошибка. Я исправил это в своем описании. Большое спасибо за ответ и профессию!
Сумма также немного подозрительна, потому что количество слоев в месяц будет варьироваться.
Извините, это последний код 'x=tapp(r,ym,sum)'?