Параллельная обработка изображений MODIS в R

Question

Я загрузил все данные о растительности из MODIS для Африки, и теперь я пытаюсь создать мозаики из изображений и сохранить их как geo tif. У меня нет проблем делать все это, пока я не попытаюсь запустить его в параллелях на MacPro 6cores, 12 потоков. Код работает, но он использует только 1% от доступности ядра и требует навсегда для завершения процессов.

Мне очень нужна помощь, у меня есть более 70 ГБ грамм 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])) 
} 

Answer 1

Зачем использовать R вообще, когда вы можете сделать это в две строки, из оболочки, с многопоточной и быстрее с помощью GDAL utilities?

Если я правильно понял ваш вопрос и код, вы хотите извлечь подтаборку NDVI из, предположительно, контейнеров MOD13Q1, HDF5, перепрограммировать на lat/lon, изменить тип файла на GeoTiff и мозаику всех гранул, охватывающих Африку, на одну большую мозаику.

1. Сборка a Virtual Raster Table, содержащая все подкатегории NDVI.

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

Обратите внимание на использование -sd 1, который принимает только первый subdataset из каждого контейнера HDF.

2. Reproject и изменить формат.

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.

Answer 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, который также позволяет использовать многопоточность и, следовательно, немного ускоряет вычисления.

Надеюсь, это полезно.