Объединение нескольких растров в R

Я пытался найти эффективный по времени способ объединить несколько растровых изображений в R. Это смежные сцены ASTER из южного региона Килиманджаро, и моя цель - собрать их вместе, чтобы получить одно большое изображение.

Вот что у меня получилось (объект ast14dmo, представляющий список объектов RasterLayer):

# Loop through single ASTER scenes
for (i in seq(ast14dmo.sd)) {
  if (i == 1) {
    # Merge current with subsequent scene
    ast14dmo.sd.mrg <- merge(ast14dmo.sd[[i]], ast14dmo.sd[[i+1]], tolerance = 1)
  } else if (i > 1 && i < length(ast14dmo.sd)) {
    tmp.mrg <- merge(ast14dmo.sd[[i]], ast14dmo.sd[[i+1]], tolerance = 1)
    ast14dmo.sd.mrg <- merge(ast14dmo.sd.mrg, tmp.mrg, tolerance = 1)
  } else {
    # Save merged image
    writeRaster(ast14dmo.sd.mrg, paste(path.mrg, "/AST14DMO_sd_", z, "m_mrg", sep = ""), format = "GTiff", overwrite = TRUE)
  }
}

Как вы наверняка догадались, код работает. Однако слияние занимает довольно много времени, учитывая, что размер каждого отдельного растрового объекта составляет около 70 МБ. Я также пробовал Reduce и do.call, но это не удалось, так как я не мог передать аргумент «толерантность», который позволяет обойти различное происхождение растровых файлов.

Кто-нибудь знает, как ускорить процесс?


person fdetsch    schedule 08.04.2013    source источник


Ответы (7)


Вы можете использовать do.call

ast14dmo.sd$tolerance <- 1
ast14dmo.sd$filename <- paste(path.mrg, "/AST14DMO_sd_", z, "m_mrg.tif", sep = "")
ast14dmo.sd$overwrite <- TRUE
mm <- do.call(merge, ast14dmo.sd)

Вот некоторые данные из примера в raster::merge

r1 <- raster(xmx=-150, ymn=60, ncols=30, nrows=30)
r1[] <- 1:ncell(r1)
r2 <- raster(xmn=-100, xmx=-50, ymx=50, ymn=30)
res(r2) <- c(xres(r1), yres(r1))
r2[] <- 1:ncell(r2)

x <- list(r1, r2)
names(x) <- c("x", "y")
x$filename <- 'test.tif'
x$overwrite <- TRUE
m <- do.call(merge, x)
person Robert Hijmans    schedule 15.04.2013
comment
Также отличное решение, спасибо! Я только что взглянул на время вычислений, и оказалось, что ваш подход через do.call работает почти в два раза быстрее, чем Reduce. - person fdetsch; 16.04.2013
comment
По какой-то причине это не работает с именованным списком, что может вызвать проблемы. - person cmbarbu; 16.08.2019
comment
В этом случае вам нужно добавить эту строку, я думаю, names(x) <- c("x", "y") - person Robert Hijmans; 16.08.2019

Функция «слияния» из пакета Raster работает немного медленнее. Для больших проектов более быстрый вариант - работать с командами gdal в R.

library(gdalUtils)
library(rgdal)

Создайте список всех растровых файлов, к которым вы хотите присоединиться (в вашем текущем рабочем каталоге).

all_my_rasts <- c('r1.tif', 'r2.tif', 'r3.tif')

Создайте растровый файл шаблона для построения. Подумайте об этом как о большом пустом холсте, на который можно добавить плитки.

e <- extent(-131, -124, 49, 53)
template <- raster(e)
projection(template) <- '+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs'
writeRaster(template, file="MyBigNastyRasty.tif", format="GTiff")

Объедините все растровые плитки в один большой растр.

mosaic_rasters(gdalfile=all_my_rasts,dst_dataset="MyBigNastyRasty.tif",of="GTiff")
gdalinfo("MyBigNastyRasty.tif")

Это должно работать очень хорошо для скорости (быстрее, чем слияние в растровом пакете), но если у вас есть тысячи тайлов, вы можете даже сначала изучить создание vrt.

person Matthew Bayly    schedule 18.09.2016
comment
Нет ли более эффективного способа установить экстент, не угадывая широту и долготу? - person Herman Toothrot; 06.01.2018
comment
@HermanToothrot: не глядя на каждый отдельный файл, чтобы узнать окончательный размер объединенных файлов, вам нужно либо сделать (щедрое) предположение (образованное, например, открыв файлы границ в ГИС). Чтобы получить точные экстенты окончательного объединенного файла, вы можете выполнить быстрое объединение экстентов, прежде чем выполнять настоящее объединение. FINAL_EXTENT = Reduce(function(d1, d2) merge(extent(raster(d1)), extent(raster(d2))), RASTER_LIST) - person Honeybear; 23.02.2021
comment
одним из недостатков mosaic_rasters, по-видимому, является то, что невозможно определить поведение перекрытия. поведение, похоже, заключается в сохранении ВСЕХ значений первого растра, даже если они являются значениями NA. - person Honeybear; 23.02.2021

Вы можете использовать Reduce, например, так:

Reduce(function(...)merge(...,tolerance=1),ast14dmo.sd)
person agstudy    schedule 08.04.2013

Инструмент мозаики SAGA GIS (http://www.saga-gis.org/saga_tool_doc/7.3.0/grid_tools_3.html) дает максимальную гибкость для объединения числовых слоев и по умолчанию выполняется параллельно! Вам нужно только сначала перевести все растры / изображения в формат SAGA .sgrd, а затем запустить командную строку saga_cmd.

person Tom Hengl    schedule 05.07.2019

Я столкнулся с той же проблемой, и я использовал

#Read desired files into R
data_name1<-'file_name1.tif' 

r1=raster(data_name1)

data_name2<-'file_name2.tif'

r2=raster(data_name2)

#Merge files
new_data <- raster::merge(r1, r2)

Хотя он не создавал новый объединенный растровый файл, он сохранялся в среде данных и создавал объединенную карту при построении.

person Sijeh Asuk    schedule 11.02.2019

Я столкнулся со следующей проблемой, когда пытался наложить мозаику из нескольких растров друг на друга.

In vv[is.na(vv)] <- getValues(x[[i]])[is.na(vv)] :
  number of items to replace is not a multiple of replacement length 

Как отметил @Robert Hijmans, это, вероятно, из-за смещения растров. Чтобы обойти это, мне пришлось сначала пересэмплировать растры.

library(raster)

x  <- raster("Base_raster.tif")
r1 <- raster("Top1_raster.tif")
r2 <- raster("Top2_raster.tif")

# Resample
x1 <- resample(r1, crop(x, r1))
x2 <- resample(r2, crop(x, r2))

# Merge rasters. Make sure to use the right order
m <- merge(merge(x1, x2), x)

# Write output
writeRaster(m,
            filename = file.path("Mosaic_raster.tif"),
            format = "GTiff",
            overwrite = TRUE)
person Tung    schedule 26.12.2020

Я протестировал решение с использованием gdalUtils, предложенного Мэтью Бейли. Работает довольно хорошо и быстро (у меня есть около 1000 изображений, которые нужно объединить). Однако после проверки документа mosaic_raster функции здесь , Я обнаружил, что он работает без создания растрового шаблона перед мозаикой изображений. Я вставил примеры кодов из документа ниже:

outdir <- tempdir()
gdal_setInstallation()
valid_install <- !is.null(getOption("gdalUtils_gdalPath"))
if(require(raster) && require(rgdal) && valid_install)
{
layer1 <- system.file("external/tahoe_lidar_bareearth.tif", package="gdalUtils")
layer2 <- system.file("external/tahoe_lidar_highesthit.tif", package="gdalUtils")
mosaic_rasters(gdalfile=c(layer1,layer2),dst_dataset=file.path(outdir,"test_mosaic.envi"),
    separate=TRUE,of="ENVI",verbose=TRUE)
gdalinfo("test_mosaic.envi")

}

person Feng Tian    schedule 02.01.2021