在 R 中并行处理 MODIS 图像答案

【问题标题】：Parallel processing MODIS images in R在 R 中并行处理 MODIS 图像
【发布时间】：2015-10-12 17:12:04
【问题描述】：

我已经从 MODIS 下载了非洲的所有植被数据，现在我正在尝试从图像中创建马赛克并将它们保存为 geo tif。在我尝试在 MacPro 6 核、12 线程上并行运行它之前，我做这一切都没有问题。代码可以运行，但它只使用了 1% 的核心可用性，并且需要永远完成这些过程。

我真的需要帮助，我需要将超过 70GB 的 MODIS 颗粒转换为马赛克，如果我不使用实验室中获得的所有计算机能力，这将永远需要。

这里是并行处理的代码部分：

有两个部分并行运行，1.选择NDVI； 2. 要创建马赛克，将它们投影并保存为 .tif。

# register the cluster with 10 cores 
registerDoParallel(cores=11)
miCluster<-makeCluster(11)
registerDoParallel(miCluster)
library(doMC)
registerDoMC(11)

# 1.
# select the NDVI for each Africa Granule and put it in SdsList
sdsList <- foreach(k = 1:length(dateGranules), .packages=c("raster", "gdalUtils","foreach"))%dopar%{
  for(j in 1:1:length(dateGranules)){
    return(sapply(X=dateGranules[[j]], FUN=function(x){get_subdatasets(x)[1]}))
  }
}

2.
# Generate the Mosaic for Africa with NDVI as aoutput
foreach(j = 1:length(sdsList), .packages=c("raster", "gdalUtils", "foreach"))%dopar%{
    gdalwarp(srcfile=sdsList[[j]], t_srs="+proj=longlat +datum=WGS84   +no_defs",
      dstfile=file.path(dest, names[j]))
}

【问题讨论】：

标签： r foreach parallel-processing gdal

【解决方案1】：

当您可以在两行中使用 R 时，为什么要使用 R，从 shell，多线程并且使用 GDAL utilities 更快？

如果我正确理解您的问题和代码，您想从大概是 MOD13Q1、HDF5 容器中提取 NDVI 子数据集，重新投影到纬度/经度，将文件类型更改为 GeoTiff 并将覆盖非洲的所有颗粒镶嵌成一个大马赛克。

构建一个包含所有 NDVI 子数据集的 Virtual Raster Table。

gdalbuildvrt -sd 1 NDVI.vrt MOD13Q1.*.hdf

注意-sd 1 的使用，它只从每个 HDF 容器中获取第一个子数据集。

重新投影并更改格式。

gdalwarp gdalwarp -t_srs "+proj=longlat +datum=WGS84 +no_defs" -multi NDVI.vrt NDVI-mosaic.tif

第一步几乎不需要任何时间，因为它只创建一个元数据文件，其中包含有关 HDF 文件层次结构的信息，并且没有进行任何实际处理。如果你的颗粒是相邻的，例如h15v05, h15v04, h14v05,...，它们将被它们的地理位置正确镶嵌。

第二步是实际的重投影和文件格式更改。我在我的 core2duo (3GHz) 上对此进行了测试，重新投影 20 个 MODIS 颗粒大约需要 1.5 分钟。所以即使没有多线程，你也应该比 R 快很多。

【讨论】：

是的，这正是我想做的。一旦我从流感中恢复过来，我就会照你说的做。非常感谢！
Kersten，我按照你的建议做了，但我收到以下错误，我拥有的任何数据集都不会超过 40% 的进度：0...10...20...30... 40.ERROR 4: `HDF4_EOS:EOS_GRID:"MOD13A3.A2005213.h20v08.005.2008069184257.2005-08-01.hdf":MOD_Grid_monthly_1km_VI:1 km 月 NDVI' 在文件系统中不存在，并且未被识别为受支持数据集名称。
您的 R 方法是否出现同样的错误？它是在构建 VRT 时发生还是在变形时发生？听起来数据集与其他 MODIS hdf 容器不同，甚至可能已损坏。
R 中没有发生该错误。我什至下载了似乎已损坏的文件。我认为这是一个 gdal 问题，因为我必须启用 hdf4 功能，因为 brew 版本只附带 hdf5。你怎么看？
很可能是 gdal 问题，尤其是 HDF。你到底使用的是什么数据集，它是在 HDF4 中提供的？据我所知，目前所有的 5 个数据集都是 HDF5。

【解决方案2】：

只是对这个问题的一点更新。

如果您仍然需要坚持使用 R 来处理您的图像，您可以使用 gdalUtils 包，它基本上是所有 GDAL 应用程序的包装器。

例如，您可以使用 @Kersten 在 R 中巧妙地提出的命令，方法是使用以下脚本：

library(gdalUtils)

gdalbuildvrt("MOD13Q1.*.hdf", "NDVI.vrt", sd=1)

gdalwarp("NDVI.vrt", "NDVI-mosaic.tif",
         t_srs="+proj=longlat +datum=WGS84 +no_defs", multi=TRUE)

这些函数基本上会调用你系统的GDAL命令，但是通过R。

请注意multi 参数，它还允许多线程，因此可以稍微加快计算速度。

希望这是有用的。

【讨论】：