hunzikp / velox

https://hunzikp.github.io/velox/
119 stars 23 forks source link

Access cell number of extracted value? #33

Open mikejohnson51 opened 5 years ago

mikejohnson51 commented 5 years ago

Hi,

Thanks for such a great package. I am looking for a way to access the cell number(or lat long) of the values extracted using velox_extract.

The goal is that the values can be modified and then replaced.

I've looked through the examples on your website and the documentation but haven't seen anything.

Ideally, it would return a similar result as ...

# pl is a polyline
# r is a raster
elevations <- raster::extract(r, pl, cellnumber = T)[[1]]

head(elevations)

   cell    layer
1 1959129 22.73608
2 1960496 28.27165
3 1961863 35.43027
4 1963230 35.97610
5 1963231 35.62851
6 1964598 36.00765

but with all the benefits of your package!

Either the cell index or lat/lon would be awesome if you know a solution.

Thank you!

Mike

jdonager commented 5 years ago

I'd also be very interested in this. Has anyone come up with a workaround? Thanks!

mikejohnson51 commented 5 years ago

Hi Jonathon,

I never got this working with velox but the following workflow is pretty quick using base raster and was enough for me to process state wide elevation and river network data in California.

library(sp)
library(raster)

# Base CRS
crs <- sp::CRS("+init=epsg:4326")

# End Points of a Line
coords <- rbind(c(-100, -90), c(80, 90))
line <- as(sp::SpatialPoints(coords, proj4string = crs),  "SpatialLines")

# Random Raster
r <- raster::raster(nrows = 100, ncols = 100, ymn = -80, ymx = 80, crs = crs)
r[] <- runif(raster::ncell(r))

# Get data.frame with intersecting cell number and value
v <- data.frame(extract(r, line, cellnumber = T, along = T)[[1]])

# Extract row and col indices
df = data.frame(rowColFromCell(r, v$cell), v$cell, v$layer)

# Rename Cols
names(df) = c('row', 'col', 'cell', 'val')

# Print
head(df)

  row col cell       val
1 100  25 9925 0.9834805
2 100  26 9926 0.2279086
3  99  26 9826 0.3707274
4  98  26 9726 0.1424308
5  98  27 9727 0.3619351
6  97  27 9627 0.9143546

Hope it helps!

Mike

ghost commented 5 years ago

I have found a solution with velox, supposing you wanted to extract values from a raster r using a spatial* object called line: `

First create a new raster object with the cell numbers

rID=r rID[]=1:(nrow(rID)*ncol(rID))

Stack the two rasters

r=stack(rID,r)

Create a velox object

r=velox(r)

Finally extract the values from the rasterstack

r$extract(sp=line) ` Hope it helps! Fabrice