Hi -
I have a 2 meter DEM, a 30m satellite image and some point data. Basically,
I would like to extract the 2m dem pixels that correspond with the 30m
satellite data, but only for the satellite pixels that correspond with a
point observation. i've been trying to get at the answer from the sp
vignette 'over' but am stuck.
an example:
#define DEM
dem2m=raster(x=matrix(runif(22500,min=2000,max=2500),nrow=150),crs='+proj=utm
+zone=13 +datum=NAD83')
extent(dem2m)=c(0,300,0,3000)
#
#define satellite image
fsca.dem=raster(x=matrix(runiff(100,0,100),nrow=10),crs='+proj=utm +zone=13
+datum=NAD83')
extent(dem2m)=c(0,300,0,3000)
#
#define point data
obs=data.frame(x=runif(5,0,300),y=runif(5,0,300),sd=runif(5,0,5))
coordinates(obs)=~x+y
proj4string(obs)='+proj=utm +zone=13 +datum=NAD83'
#
#--- in my case, I resampled the sat image to the dem grid. but thats
already taken care of in the example.
##dem.pe=projectExtent(dem2m,crs=projection(dem2m))
##res(dem.pe) <- 30
##fsca.dem=projectRaster(fsca,dem.pe,method='ngb')
#
# This is where I am lost
# ----------- I thought I'd use polygons to extract the dem data. The catch
is, rather than extracting every polygon of data (which takes longer than
I've been willing to wait given 2e6 pixels) I would like to subset this
spatialpolygondf to the 5 polygons that have observational data. if the
same polygon has multiple obs in it, it should be replicated.
fsca.poly=rasterToPolygons(fsca.dem)
#first i tried this but this extracts the fsca values. I need this too
eventually but was hoping for a subsetted spatialpolygon.dataframe
fsca.pts=over(ss.depth,fsca.poly)
#second i tried this but the list has rownames that I can't interpret.
(fsca cellnumbers?)
fsca.pts.poly=over(ss.depth,fsca.poly,returnList=T)
polynames=unlist(lapply(fsca.pts.poly,rownames))
#but what i'd like to do is apply the cv function on the 225 2m elev pixels
that correspond to the 30 fsca pixel for which i have an observation point.
elev.cv=extract(dem2m,*y=fsca.poly.w.pts*,fun=cv)
# as a note, the buffer option in extract would not work currently because
the point data isn't necessarily in the middle dem pixel relative to the
fsca pixel. so one solution might be to 1. subset the polygon dataframe 2.
use the cell center from the fsca pixel to extract from the dem raster with
a buffer of 7 cells.
anyway, i appreciate any help you can give, in particulary with subsetting
the spatialpolygondf.
Thanks
Dominik
R version 3.0.2 (2013-09-25)
Platform: x86_64-apple-darwin10.8.0 (64-bit)
locale:
[1] C
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] maptools_0.8-30 reshape2_1.2.2 plyr_1.8
rgdal_0.8-16
[5] ggplot2_0.9.3.1 RColorBrewer_1.0-5 scales_0.2.3
raster_2.2-31
[9] sp_1.0-15
loaded via a namespace (and not attached):
[1] MASS_7.3-29 colorspace_1.2-4 compiler_3.0.2 dichromat_2.0-0
[5] digest_0.6.4 foreign_0.8-55 grid_3.0.2 gtable_0.1.2
[9] labeling_0.2 lattice_0.20-29 munsell_0.4.2 proto_0.3-10
[13] stringr_0.6.2 tcltk_3.0.2 tools_3.0.2
[[alternative HTML version deleted]]