CODEWITHSUNDEEP
rasterr-rasterrasterizingrasterize
ctlamb
ctlamb

Reputation: 133

Rasterize() slow on large SpatialPolygonsDataFrame, alternatives?

I have a large (266,000 elements, 1.7Gb) SpatialPolygonsDataFrame that I am try to convert into 90m resolution RasterLayer (~100,000,000 cells)

The SpatialPolygonsDataFrame has 12 variables of interest to me, thus I intend to make 12 RasterLayers

At the moment, using rasterize(), each conversion takes ~2 days. So nearly a month expected for total processing time.

Can anyone suggest a faster process? I think this would be ~10-40x faster in ArcMap, but I want to do it in R to keep things consistent, and it's a fun challenge!

general code

######################################################
### Make Rasters
######################################################
##Make template
r<-raster(res=90,extent(polys_final))

##set up loop
loop_name <- colnames(as.data.frame(polys_final))

for(i in 1:length(loop_name)){
a <-rasterize(polys_final, r, field=i)
writeRaster(a, filename=paste("/Users/PhD_Soils_raster_90m/",loop_name[i],".tif",sep=""), format="GTiff")
}

Upvotes: 2

Views: 1711

Answers (1)

Andrew Barr
Andrew Barr

Reputation: 3800

I think this is a case for using GDAL, specifically the gdal_rasterize function.

You probably already have GDAL installed on your machine if you are doing a lot of spatial stuff, and you can run GDAL commands from within R using the system() command. I didn't do any tests or anything but this should be be MUCH faster than using the raster package in R.

For example, the code below creates a raster from a shapefile of rivers. This code creates an output file with a 1 value wherever a feature exists, and a 0 where no feature exists. 

path_2_gdal_function <- "/Library/Frameworks/GDAL.framework/Programs/gdal_rasterize"
outRaster <- "/Users/me/Desktop/rasterized.tiff"
inVector <- "/Full/Path/To/file.shp"
theCommand <- sprintf("%s -burn 1 -a_nodata 0 -ts 1000 1000 %s  %s", path_2_gdal_function, inVector, outRaster)
system(theCommand)
  • the -ts argument provides the size of the output raster in pixels
  • the -burn argument specifies what value to put in the output raster where the features exist
  • -a_nodata indicates which value to put where no features are found

For your case, you will want to add in the -a attribute_name argument, which specifies the name of the attribute in the input vector to be burned into the output raster. Full details on possible arguments here.

Note: that the sprintf() function is just used to format the text string that is passted to the command line using the system() function

Upvotes: 3

Related Questions