2

I am attempting to use R to rasterize some downloaded shapefiles that include polygons with holes in them and R is failing to recognize the holes.

The holes seem to be designated as hole=T, so I don't think the issue is the one solved in the question: rasterize ESRI shapefile with holes but FALSE hole slots

What else could be going on? How do I troubleshoot this issue?

I've been using the rasterize function without any changes to the default settings, and the data I've been working with are reprojected from shapefiles here: ftp://ftp.granit.sr.unh.edu/pub/GRANIT_Data/Vector_Data/Environment_and_Conservation/d-wap/, such as hemhwdpine.shp.

Community
  • 1
  • 1
user1521655
  • 486
  • 3
  • 17
  • Have you figured out which of the 20598 polygons in that file are causing you problems? Doing so would help us help you (if you don't happen to figure out the issue while narrowing down the problem yourself). – Josh O'Brien Apr 30 '15 at 00:24
  • The polygon with FGID 19334 is one offender, but I'm certain that it isn't the only one. – user1521655 Apr 30 '15 at 02:25
  • @Josh O'Brien, I just read your response to this question: http://stackoverflow.com/questions/14992197/shapefile-to-raster-conversion-in-r?rq=1, and realized that I have read in my shapefiles using maptools rather than rgdal. Could that be part of the problem? Why is rgdal the better choice? – user1521655 Apr 30 '15 at 02:43
  • One reason `rgdal::readOGR()` is better is that it reads in and retains the layer's CRS information. To see that for yourself, compare the value returned by `proj4string()` applied to objects read in by `maptools::readShapePoly()` and `rgdal::readOGR()`. I actually read the file in using `raster::shapefile()`, which is a convenient wrapper around `readOGR()`, and had no trouble rasterizing the complicated polygon with FGID==19334, so it looks like using `readOGR()` may be enough to solve your problem. – Josh O'Brien Apr 30 '15 at 06:24
  • Tried switching and still not getting the holes. Wondering now if it might have to do with limited disk space for temporary files, or maybe a memory problem of some sort. – user1521655 Apr 30 '15 at 15:14
  • Is the resolution of your raster perhaps too coarse to pick up the holes? I had the same initial thought about memory limitations, but one of the nice things about **raster** as a package is that it does not depend on holding really large rasters in memory. Even if you do not supply `rasterize()` with a `filename` argument, it checks internally whether `canProcessInMemory(rstr,3)` and if not, sets up a tempfile to work off of. (See `raster:::.polygonsToRaster` for the exact code that's called). – Josh O'Brien Apr 30 '15 at 15:42
  • That said, I've had problems with `raster::rasterize()` being way too slow when producing huge rasters, to such an extent that I've sometimes turned to `gdalUtils::gdal_rasterize()`. It's a very thin wrapper around [`gdal_rasterize`](http://www.gdal.org/gdal_rasterize.html) and, though it requires a little bit more setup on my part, it can be much much faster. – Josh O'Brien Apr 30 '15 at 15:45
  • Let us [continue this discussion in chat](http://chat.stackoverflow.com/rooms/76687/discussion-between-user1521655-and-josh-obrien). – user1521655 May 01 '15 at 00:51
  • Update: the commands worked fine on a different computer. Problem not diagnosed, but clearly has something to do with my system/configuration. – user1521655 May 01 '15 at 13:58
  • Aah, that's good news. Thanks for following up. – Josh O'Brien May 01 '15 at 14:11

0 Answers0