0

I'm trying to convert a polygon shapefile with small values. Values from the colunm propEmp range from 0.000002 to 0.119419. This is my attempt:

# Load shapefile
emp_planejado <- shapefile("./planejado/7_wgs/emp_planejado.shp")

# Load raster model
r_bioma<- raster("./_GRID/grid_caatinga_disol_64bit.tif")

# List values from tipologia field
tipologia<-unique(emp_planejado$tipologia)


for (tp in tipologia){
# Select features for each value in tipologia
tipo<- emp_planejado[emp_planejado$tipologia==tp,]

# Rasterize
r_pol <- rasterize(tipo,r_bioma,field="propEmp",background=NA,mask=F)

# Merge 
raster_merge <- merge(r_pol,r_bioma)

# Save raster
writeRaster(raster_merge,filename= paste0("./planejado/8_raster/",tp,"_planejado"),format="GTiff",NAflag=-9999,overwrite=TRUE)  
}

r_bioma is a 64bit with double precision raster with all pixel values equal to 0.

There is no overlap between features of this polygon, just boundary touching, so I did not use funas an argument of rasterize.

After rasterize, when I check the minValue and maxValue from r_pol, instead of getting 0.000002 and 0.119419, I get 0.08687903 and 0.1140689.

I don't know where the problem is. Could you help me?

Tiago
  • 143
  • 1
  • 2
  • 7
  • You could try with `gdalUtils::gdal_rasterize`, e.g. `gdal_rasterize('planejado/7_wgs/emp_planejado.shp', '_GRID/grid_caatinga_disol_64bit.tif', a='propEmp')`. – jbaums Sep 18 '15 at 23:15

1 Answers1

1

Your question is not clear, your example not reproducible, and we do not know what you are trying to achieve. So it is very difficult to provide meaningful help.

First, r_pol is created inside of a loop, in which tipo is subsetted, so that it would be expected that the range of values are not the same.

More fundamentally, why are you doing this complex loop? It would seem that what you want is something like the below, but I cannot know that as you did not say what you wanted to achieve.

library(raster)
emp_planejado <- shapefile("./planejado/7_wgs/emp_planejado.shp")
r_bioma <- raster("./_GRID/grid_caatinga_disol_64bit.tif")

r_pol <- rasterize(emp_planejado, r_bioma, field="propEmp", background=0, filename="./planejado/8_raster/planejado.tif",overwrite=TRUE)

Perhaps this does not do it, but it hard to imagine that you actually would need the complex loop that you are using.

Finally, it may be that some polygons are very small and do not cover a cell. In that case, their values could get lost.

Robert Hijmans
  • 40,301
  • 4
  • 55
  • 63
  • Sorry if I've not make myself clear. The thing is: I have a polygon shapefile wtih many features and I want to make a raster layer of each one of them, that's why I did the loop. I've discovered what the problem was, and it was exactly what you have suspect: there were some polygons with a very small area, so they were ignored in rasterizing. To solve this problem, I converted the polygons to points and then rasterized the points shapefile and it worked well. All valeus,even the small ones, were taken to the right pixel perfectly. Thanks for your help! – Tiago Sep 21 '15 at 19:40