R語言遍歷資料夾求取其中所有柵格檔案的平均值

疯狂学习GIS發表於2024-06-28
R語言

  本文介紹基於R語言中的raster包，遍歷讀取多個資料夾下的多張柵格遙感影像，分別批次對每一個資料夾中的多個柵格影像計算平均值，並將所得各個結果柵格分別加以儲存的方法。

  其中，本文是用R語言來進行操作的；如果希望基於Python語言實現類似的平均值求取操作，大家可以參考Python ArcPy批次計算多時相遙感影像的各像元平均值與Python忽略NoData計算多張遙感影像的像元平均值：whitebox庫這兩篇文章。

  首先，來看一下本文所需實現的需求。如下圖所示，現有多個資料夾，其中每一個資料夾內部都含有大量的柵格遙感影像。

image

  其中，上圖中的每一個資料夾的命名都是以遙感影像的分幅條帶號為依據的。例如，開啟第一個名為47RMN的資料夾，其中均為條帶號為47RMN（即同一空間範圍）、不同成像時間的遙感影像，如下圖所示；其中，紫色框內的遙感影像檔名即可看出，這些影像是同一條帶號、不同時間的遙感影像資料。

  我們要做的，就是分別對每一個資料夾中的全部遙感影像計算平均值，從而得到不同條帶號遙感影像的平均值；最終我們將得到多張結果影像，每一景結果影像就是這一條帶號、不同成像時間對應的遙感影像的平均值。同時為了方便區分，我們需要將每一景結果影像檔案的檔名設定為與條帶號有關的內容。

  明確了需求，我們即可開始程式碼的撰寫。本文所用到的程式碼如下所示。

library(raster)
result_path <- r"(E:\02_Project\01_Chlorophyll\Select\Result)"
tif_folder <- list.files(path = r"(E:\02_Project\01_Chlorophyll\Select)", pattern = NULL, all.files = FALSE, full.names = TRUE)
for (folder in tif_folder){
  folder_name <- substr(folder, nchar(folder) - 4, nchar(folder))
  tif_file_name <- list.files(path = folder, pattern = ".tif$", full.names = TRUE, ignore.case = TRUE)
  tif_file_all <- stack(tif_file_name)
  NAvalue(tif_file_all) <- -10000
  tif_mean <- calc(tif_file_all, fun = mean, na.rm = TRUE)
  tif_mean_new <- tif_mean / 100
  # plot(tif_mean_new)
  result_file_name <- file.path(result_path, paste(folder_name, "_mean.tif", sep = ""))
  rf <- writeRaster(tif_mean_new, filename = result_file_name, overwrite = TRUE)
  cat(folder_name, "is completed!", "\n")
}

  首先，需要透過library(raster)程式碼，匯入本文所需的R語言raster包；關於這一包的配置，大家可以參考基於R語言的raster包讀取遙感影像。接下來，我們需要指定結果存放的路徑，並將其放入變數result_path中。

  接下來，我們透過list.files()函式，將包含有各個條帶號的小資料夾的大資料夾（也就是本文開頭第一張圖所示的資料夾）加以遍歷，將每一個小資料夾的路徑存入tif_folder。執行上述前3行程式碼後，得到的tif_folder結果如下圖所示。

  可以看到，tif_folder是一個字串，其中每一個元素都是每一個小資料夾的路徑。

  接下來的for迴圈，就是對tif_folder加以遍歷，即對每一個小資料夾進行操作。其中，我們首先透過substr()函式，獲取當前操作的小資料夾名稱，並將其存放於folder_name中；隨後，對當前對應的小資料夾加以遍歷，取出其中的全部遙感影像檔案，並存放於tif_file_name；接下來，就是讀取全部遙感影像，並計算其平均值；這裡具體的程式碼解釋大家可以參考文章R語言求取大量遙感影像的平均值、標準差：raster庫。此外需要注意的是，由於我這裡每一景遙感影像原本沒有專門設定NoData值，而是用-10000作為其NoData值，因此需要透過NAvalue(tif_file_all) <- -10000這句程式碼，將值為-10000的像元作為NoData值的像元，防止後期計算平均值時對結果加以干擾。

  接下來，我們透過file.path()函式配置一下輸出結果的路徑——其中，結果遙感影像檔案的名稱就可以直接以其所對應的條帶號來設定，並在條帶號後新增一個_mean字尾，表明這個是平均值的結果影像；但此外，這個僅僅是檔案的名字，還需要將檔名與路徑拼接在一起，才可以成為完整的儲存路徑，因此需要用到file.path()函式。最後，將結果影像透過writeRaster()函式加以儲存即可，這句程式碼的解釋大家同樣參考R語言求取大量遙感影像的平均值、標準差：raster庫這篇文章即可。

  最後，由於我們要處理的資料夾比較多，因此可以透過cat()函式輸出一下當前程式碼的執行進度。

  執行上述程式碼，我們將在指定的結果儲存路徑中看到每一個條帶號對應的平均值結果影像，如下圖所示。

  至此，大功告成。

相關文章