[R-sig-Geo] Intersection of polygons with raster

[Kent Johnson]

Thank you Michael, that is exactly what I needed. It is fairly slow for my
case but I only need to do it once :-)

And thank you for providing an example, I'm not too familiar with rasters
and would have had some trouble creating it.

Kent

On Fri, May 18, 2018 at 11:41 PM, Michael Sumner <mdsumner at gmail.com> wrote:

> Hi Kent, this is pretty straightforward with raster::extract. If speed is
> an issue you can burn the polygon id into the grid and then group values by
> that.
>
> library(raster)
> r <- raster(volcano)
>
> ## very simplistic polygon layer example
> r2 <- raster(r)
> res(r2) <- res(r2) * 20
> p <- sf::st_as_sf(as(r2, "SpatialPolygonsDataFrame"))
>
> ## a dummy threshold value
> onefoot <- 150
> ## grab the first column from extract (it gives multiple columns for
> multi-layer rasters)
> p$onefoot <- extract(r, p, fun = function(x, na.rm = TRUE) any(x >
> onefoot))[,1]
>
> plot(sf::st_geometry(p), col = p$onefoot + 1)
> contour(r, levels = onefoot, col = "white", add = T)
>
>
> ## if speed is an issue, use fasterize for a more abstract workflow
> library(fasterize)
> p$rownum <- 1:nrow(p)
> idgrid <- fasterize(p, r, field = "rownum")
> p$onefoot <- tapply(values(r), values(idgrid), function(x, na.rm = TRUE)
> any(x > onefoot))
> plot(p)
>
> FWIW it's always useful to provide a reprex, since that does take time and
> might miss the mark for your situation. It's also trickier to provide two
> objects that are relevant to each other, so any upfront work you can do in
> an example helps the answerer.
>
> HTH
>
> On Sat, 19 May 2018 at 12:20 Kent Johnson <kent3737 at gmail.com> wrote:
>
>> Hi,
>>
>> I have a raster object `flood` containing projected flooding levels and a
>> simple features object `parcels` containing MULTIPOLYGONs representing
>> property parcels. I would like to find all the parcel polygons for which
>> there is any flood > 1 foot. What function(s) can do this? Something like
>> raster::intersect(parcels, flood >= 1)
>>
>> I'll put together a reprex if that helps. Hoping someone can point me to
>> the right function to do this, I haven't found anything promising.
>>
>> Thanks,
>> Kent
>>
>>         [[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
>>
> --
> Dr. Michael Sumner
> Software and Database Engineer
> Australian Antarctic Division
> 203 Channel Highway
> <https://maps.google.com/?q=203+Channel+Highway+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
> Kingston Tasmania 7050 Australia
> <https://maps.google.com/?q=203+Channel+Highway+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
>
>

	[[alternative HTML version deleted]]