I've got a shapefile that I would like to convert into a raster, but the trick is that the attribute I want a raster version of is character data. I'm familiar with geoprocessing using GDAL (in the command line) and R, so solutions using either of those two are preferred. However, I have access to QGIS if necessary, and even ArcMap.
I'm using this shapefile. It contains polygons of municipalities in part of Brazil. The attribute NAME_2 is the name of the municipality, and I want a raster version of that.
I know that you can't make a raster with characters, so I created a numeric version of NAME_2 in R by doing
library("rgdal") munisCrop.shape = readOGR(dsn=getwd(), layer="munisCrop") # Read shapefile munisCrop.shape$NAME_2_NUM = as.numeric(munisCrop.shape$NAME_2) # Make new attribute writeOGR(obj = munisCrop.shape, dsn = getwd(), layer = "munisCrop", driver = "ESRI Shapefile") # Save new version of shapefile
munisCrop.shape$NAME_2_NUM range from 1 to 787.
I was hoping to then use
gdal_rasterize on my command line (GDAL version 1.10.1 on Mac OS X 10.9.2) to make a raster version of the shapefile. Here's the code I was trying (separate lines for clarity):
gdal_rasterize -of GTiff -ot UInt32 -a NAME_2_NUM \ -tr 0.0002777778 0.0002777778 -ot BYTE \ -co COMPRESS=LZW -l munisCrop \ munisCrop.shp munisCrop.tif
This completes easily enough. However, when viewing munisCrop.tif in QGIS, I see that the data range from 0 to 255. This is confirmed in R with
library("raster") munisCrop.raster = raster("munisCrop.tif") sort(unique(getValues(munisCrop.raster)))
0 makes sense, as NA values got converted to that. However, the maximum value of 255 is no good; remember, the maximum value of that attribute in the shapefile is 787.
I think I understand why this is happening, but I'd appreciate if somebody could really explain it to me and suggest a workaround.
I've tried various options in
gdal_rasterize instead of
-ot UInt32 (for example,
BYTE), and the result is always the same. I've also tried rasterizing it in R:
munisCrop.blankraster = raster(ext = extent(munisCrop.shp), resolution = 0.0002777778, crs = projection(munisCrop.shp) ) munisCrop.Rraster = rasterize(x = munisCrop.shp, y = munisCrop.blankraster, field = "NAME_2_NUM" )
But that seems like it takes forever.
gdal_rasterize is way faster (a matter of seconds instead of... indefinite), so I'd love a solution using that.
Thanks in advance!