CODEWITHSUNDEEP
rr-rastergeotiff
Bing-Hong Huang
Bing-Hong Huang

Reputation: 135

Using mosaic in r to merge multiple geotiff

I have 50 geotiff files in the same folder. All of them represent elevation data in different parts of the world. I would like to merge certain geotiff files, and I found mosaic in R might help us. I have moved those geotiff into the same folder, and I wrote a R script show below:

setwd()
a<-list.files(pattern="*.tiff",file.name=TRUE)
combind<-merge(a,fun=mean)

However, this script returned an error: error in as.data.frame(y)

May I ask how could I improve my script?

Upvotes: 8

Views: 8208

Answers (1)

loki
loki

Reputation: 10350

Edit 2023-06: use terra::vrt() or terra::mosaic()

Since gdalUtils seems unmaintained for a while now, you can use a similar approach as in my original answer just using the more up-to-date terra package which is the successor to the raster package:

library(terra)

vrt(
  x = list.files(path = "folder/to/images", pattern = "*.tif$", full.names = TRUE), 
  filename = "dem.vrt"
)
# afterwards read it as if it was a normal raster:
dem <- rast("dem.vrt")

Another possibility is to use terra::mosaic(). This gives you more control over the aggregation function in overlapping pixels. For that you first create a raster collection using sprc() from a list of rast objects and then mosaic to combine them into one raster.

library(terra)
# vector of file names
fls <- list.files("your/folder", ".tif$", full.names = TRUE)
# list of rast objects
r_lst <- lapply(fls, rast)
# create spatial raster collection
coll <- sprc(r_list)

# combine all rasters
mosaic(coll, function = "mean")

Original answer using gdalUtils

You can make use of the powerful GDAL functions. From my experience these are much faster than pure R code.

My approach would be with library(gdalUtils):

First, build a virtual raster file (vrt):

library(gdalUtils)
setwd(...)
gdalbuildvrt(gdalfile = "*.tif", # uses all tiffs in the current folder
             output.vrt = "dem.vrt")

Then, copy the virtual raster to a actual physical file:

gdal_translate(src_dataset = "dem.vrt", 
               dst_dataset = "dem.tif", 
               output_Raster = TRUE # returns the raster as Raster*Object
                                    # if TRUE, you should consider to assign 
                                    # the whole function to an object like dem <- gddal_tr..
               options = c("BIGTIFF=YES", "COMPRESSION=LZW"))

Another pure (and probably slower) raster package solution would be:

f <- list.files(path = "your/path", pattern = ".tif$", full.names = TRUE)
rl <- lapply(f, raster)

do.call(merge, c(rl, tolerance = 1))

you have to adjust the tolerance since the raster files will probably not have the same origin.

Upvotes: 16

Related Questions