Mosaic, or merge, rasters in R

To do a mosaic of rasters in GIS software, whether Arc or QGIS, is simple enough command-wise. However, I was running into issues about a month ago getting my rasters to mosaic properly. Something always seemed to come up wrong, so I thought I would do it in R to see how my results would differ. I’m glad I did, and so I thought I’d put this up in case anyone else out there is having the same trouble I did. Just as a note, I think the likelihood of an error increases depending on the size of your raster. Since I was using 1/3 arcsecond, i.e. 10 m spatial resolution, the rasters were fairly large.

All of my DEM data was downloaded from Earthexplorer. This mosaic would eventually become the basis for my watershed delineation, as seen in my portfolio.

Now onto the code.

I used the libraries rgdal, gdalutils, and raster. If they’re not present, use the install.packages() function. By the way, I’ve just discovered the package, easypackages. This package will allow you to load multiple packages with just one line. I recommend it.

library(easypackages)
libraries("rgdal", "gdalUtils", "raster")

So in the next line of code, I’m loading my rasters into R. All of my rasters, in this case I’m working with geotiffs, are all in the same folder for ease of loading. I’ve set my working directory to this folder.

setwd("D:/Rasters")
a <- c('1.tif', '2.tif', '3.tif', '4.tif')

The following is perhaps the most annoying part of mosaicing rasters in R. You need to create a template raster. This is basically the blank canvas that R will use to piece all of your rasters together. Think of it like a blank piece of paper. You’ll need the coordinates for the four corners of the blank piece of paper. The coordinates will be entered in the following format: x-minimum, x-maximum, y-minimum, y-maximum. You’ll use the extent() function to set the extent of your template raster.

e <- extent(-85, -83, 39, 41)
template <- raster(e)

As with all things geospatial, you need to set your projection. The projection for this particular raster dataset is NAD83, or espg 4269. If you’ve read my other tutorials, you know that you’ll be using the proj4string() and CRS() functions and that epsg information can be obtained here.

proj4string(template) <- CRS('+init=epsg:4269')

writeRaster(template, file="MiamiWatershed.tif", format="GTiff")
mosaic_rasters(gdalfile=a,dst_dataset="MiamiWatershed.tif",of="GTiff")
gdalinfo("MiamiWatershed.tif")

And you’re all set. Your mosaic raster should now be available in your working directory. Honestly, this will probably be my go to for mosaicing rasters. It seems to be faster than using Arc or Q.

setwd("D:/Raster")
a <- c('1.tif', '2.tif','3.tif','4.tif')
e <- extent(-85, -83, 39, 41)
template <- raster(e)
proj4string(template) <- CRS("+init=epsg:4269")
writeRaster(template, file="MiamiWatershed.tif", format="GTiff")
mosaic_rasters(gdalfile=a,dst_dataset="MiamiWatershed.tif",of="GTiff")
gdalinfo("MiamiWatershed.tif")

The featured image is Mosaic by Jenny Brown. It is used under a Creative Commons 2.0 License.

8 thoughts on “Mosaic, or merge, rasters in R

    1. This is a generic error that could be applicable to a wide range of problems. If you could provide the code that you used and where the error occurred, I would have a better idea of what went wrong.

  1. Hi I get this error. My guess is that list of rasters is wrong? the “a”
    mosaic_rasters(gdalfile=a,dst_dataset=”BigRaster.tif”,of=”GTiff”)
    Error in FUN(X[[i]], …) : invalid ‘file’ argument

    This is my code. naip_ndvi and naip_evi are my rasters. Img is a polygon I used to crop the original raster file before creating naip_ndvi and naip_evi.

    a <- c(naip_ndvi, naip_evi)
    img
    extent(img)

    e <- extent(img)
    template <- raster(e)
    proj4string(template) <- CRS('+proj=utm +zone=54 +south +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0')

    writeRaster(template, file="BigRaster.tif", format="GTiff", overwrite=TRUE)

    mosaic_rasters(gdalfile=a,dst_dataset="BigRaster.tif",of="GTiff")

    I hope you can help me!
    cheers and thanks for your code.

    1. The extent needs to be set to the rasters you’re trying to mosaic. You can mask the raster after you’ve mosaic’d them.

  2. Thanks for this, very helpful. When running the following code I get the below error.

    mosaic_rasters(gdalfile=Viet_NDVI,dst_dataset= “Vietnam_NDVI.tif”, of=”GTiff”)

    “Error in mosaic_rasters(gdalfile = Viet_NDVI, dst_dataset = “Vietnam_NDVI.tif”, : Some of the input files are missing: ndvi32-00.tif Some of the input files are missing: ndvi00-00.tif”

    I’m pretty sure I followed everything exactly, please let me know how to fix this.

  3. Its hard to tell whats wrong without seeing all the code. You’d have to provide all of the code. Its likely your mistake occurred before this line. Its just that it didn’t become apparent until you tried to run this function.

Leave a Reply to tK Cancel reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Google photo

You are commenting using your Google account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s