1

I would like to use R to:

  1. Read in an MrSID raster (.sid format)
  2. Clip the extent using a shapefile
  3. Save the new clipped raster as a geotiff

I am fine with steps #2 and 3, but I can't find any information about reading an MrSID file to R. Here is an example of what my workflow would look like if the raster was a .tif instead of a .sid:

# load packages
require(raster)
require(rgdal)

# Step 1: read in raster from .tif file
r <- raster("C:/.../raster.tif")
# Q: how do you do this with a .sid file?

# read in shapefile used for clip
shp <- readOGR(dsn="C:/.../folder", layer="shapefile")

# Step 2: clip raster to shapefile
r.clip <- crop(r, shp)

# Step 3: save clipped raster as tiff
writeRaster(r.clip, "C:/.../raster_clip.tif", format="GTiff")

I know this could be done relatively easily in something like QGIS, but I have a large number of shapefiles I would like to use as clipping extents and then do some follow-up data analysis, so I would like to use R to automate it.

I've found other questions regarding GDAL and MrSID (e.g. here), but none related to R, and I'm kind of a GDAL newbie so I haven't been able to successfully translate the knowledge.

  • I ended up just making a script to clip them all using QGIS, but if there's a way to do this in R in the future, I'd love to know it. – Sam Zipper Mar 30 '16 at 2:53
3

To read MrSID-files, you need the GDAL MrSID-plugin. Check if you have the MrSID-plugin correctly installed by the command gdalinfo --format MrSID in your commandprompt/shell (otherwise check here).

Make sure the rgdal-package installs the GDAL with the MrSID-plugin (you might have several versions installed). Linking the rgdal package to the correct GDAL is not easy - I'm still figuring out how to do this on my system - but see here and here for some guidelines.

If this works, use the readGDAL function {rgdal package} to read your data.

library(rgdal)
r <- readGDAL("C:/.../raster.sid")

Alternatively, you can use the gdal_translate function in the {gdalUtils package}. The advantage is that you can choose the GDAL version easily with gdal_setInstallation.

library(gdalUtils)
gdal_setInstallation(search_path =...) # choose GDAL version with this function
r <- gdal_translate(src_dataset = "C:/.../raster.sid",
           dst_dataset = "raster.tif",
           output_Raster = TRUE) # create R raster object

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.