[R-sig-Geo] over() in sp package

[Edzer Pebesma]

On 04/14/2015 01:27 AM, Kenny Bell wrote:
> Hi,
> 
> When I call:
> 
> over(spPolyDf, spGridDf)
> 
> The function seems to be returning the spGridDf$band1 value for the first
> grid cell that is contained in each polygon.

This is documented; you can use returnList = TRUE to get the list of
values for each polygon. Since you want to aggregate these anyway, you
might as well use

aggregate(spGridDF, spPolyDf, mean)

this is gives a simple average of pixels whose centers fall inside the
polygons. Do you need a weighted average, weighted by how much of the
pixel overlaps with each polygon?

> 
> What I would like is an overlap-area weighted average (or some provided
> function like sum) of all the grid cells contained (or partially contained)
> in each polygon element of spPolyDf.
> 
> I have code that does this by eventually using joinPolys() from PBSMapping,
> but it's slow and would be much cleaner to have a sp-only solution. Is
> there any way to do this? Some example code is below:
> 
> library(sp)
> library(rgdal)
> 
> # Download a USA shapefile
> tmpFile <- tempfile()
> download.file("
> http://www2.census.gov/geo/tiger/GENZ2013/cb_2013_us_state_20m.zip",
> destfile = tmpFile)
> unzipped <- unzip(zipfile = tmpFile, exdir = getwd())
> 
> # This is the actual reading stage
> states <- readOGR(getwd(), layer="cb_2013_us_state_20m")
> 
> # Remove all the files on disk
> file.remove(tmpFile)
> file.remove(unzipped)
> 
> tmpFile <- tempfile()
> 
> url <- paste0("http://www.prism.oregonstate.edu/fetchData.php?type",
>               "=bil&kind=recent&elem=tmax&range=daily&temporal=",
>               "20140416")
> 
> # Download the file to fileName
> download.file(url = url, mode = "wb", destfile = tmpFile)
> 
> # unzip the file to the working directory.
> unzippedFiles <- unzip(zipfile = tmpFile, overwrite = TRUE)
> 
> # read into R
> tMaxGrid <- readGDAL(fname = unzippedFiles[1])
> 
> # delete the file
> file.remove(tmpFile)
> file.remove(unzippedFiles)
> 
> # This doesn't work correctly
> states$tmax20140416 <- over(states, tMaxGrid)

this would assign a data.frame to what you want to be a vector (in a
data.frame), so it would have to be:

states$tmax20140416 <- over(states, tMaxGrid)[[1]]

I'll try to modify sp so that it'd warn you in this case.


-- 
Edzer Pebesma
Institute for Geoinformatics (ifgi),  University of Münster,
Heisenbergstraße 2, 48149 Münster, Germany; +49 251 83 33081
Journal of Statistical Software:   http://www.jstatsoft.org/
Computers & Geosciences:   http://elsevier.com/locate/cageo/
Spatial Statistics Society http://www.spatialstatistics.info

-------------- next part --------------
A non-text attachment was scrubbed...
Name: signature.asc
Type: application/pgp-signature
Size: 490 bytes
Desc: OpenPGP digital signature
URL: <https://stat.ethz.ch/pipermail/r-sig-geo/attachments/20150414/4b4b3983/attachment.bin>