[R-sig-Geo] R version of gdal_polygonize.py

[Matteo Mattiuzzi]

Hi Jonathan,
for what I understand from your question the main issue is the resizing of the pixels to a larger size. The solution below is 1 not very performant and 2 not fully memory save as polygonToRaster is used again, but maybe you find something usefull. 
Matteo


# largest area function (is definitively slow!)
maj <- function(x) {as.numeric(names(which.max(tapply(INDEX=x[,1],X=x[,2],FUN=sum,na.rm=TRUE))[1]))} # important, this function returns NULL if in the extraction area there is no raster value available!


# source raster
r   <- raster(nrow=20, ncol=20, xmn=-15, xmx=17, ymn=40, ymx=50)
r[] <- round(runif(ncell(r),1,2))


# target raster
trs <- "+proj=sinu +lon_0=0 +x_0=0 +y_0=0 +a=6371007.181 +b=6371007.181 +units=m +no_defs"
t   <- projectExtent(r,trs)
dim(t) <- c(10,10,1)
poly <- rasterToPolygons(t) # here maybe use the gdal_polygonize.py instead


# beginCluster()
# mj <- extract(r,poly,fun=maj,weights=T) # would be nice if it would work like that, but 'fun' argument must be a vector and with weights=T it becomes a matrix, I think thats the reason! So:
mj  <- extract(r, poly,weights=T) # You can make this also memory save by extracting only a sequence in poly
t[] <- sapply(mj, FUN=maj)


writeRaster(r,"./r.tif",NAflag=0)
writeRaster(t,"./t.tif",NAflag=0)
shapefile(poly,filename="./pols.shp")




>>> Jonathan Greenberg  02/12/13 3:41 PM >>>
R-sig-geo'ers:

Continuing on my quest to *precisely* reproject/resize a categorical raster
where the output cell is the category which has the maximum areal cover of
that cell (this is not the same as nearest neighbor), it is looking more
and more like the solution is to first convert a class raster to a polygon
surface, reproject/resize, then re-rasterize the output.

As such, is there a function within R that does this?  I see raster has the
"rasterToPolygons" but as far as I can tell, this produces one polygon per
cell (e.g. it doesn't merge the cells if they contain the same values --
Robert, correct me if I'm wrong here) -- this could cause some memory
issues, I imagine, on some big files.  I'm looking for something like
what gdal_polygonize.py
does.

Any thoughts?

--j

-- 
Jonathan A. Greenberg, PhD
Assistant Professor
Global Environmental Analysis and Remote Sensing (GEARS) Laboratory
Department of Geography and Geographic Information Science
University of Illinois at Urbana-Champaign
607 South Mathews Avenue, MC 150
Urbana, IL 61801
Phone: 217-300-1924
http://www.geog.illinois.edu/~jgrn/
AIM: jgrn307, MSN: jgrn307 at hotmail.com, Gchat: jgrn307, Skype: jgrn3007

    [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
R-sig-Geo at r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo