2016-02-12 87 views
3

我有一個包含在單個stack中的柵格圖層以及包含點座標的SpatialPoints對象。如果單元格包含一個點,我試圖將柵格圖層的值更改爲NA。如果單元格包含點,則將光柵單元格值更改爲NA

我在下面提供了一個可重現的例子。然而,爲了說明這一點,我的真實數據由10層「棲息地」(海拔,樹冠覆蓋等)組成,這些層都包含在堆棧中。 rasters覆蓋面積約爲34 x 26公里。另外,我有大約10,000個來自動物的GPS位置。

因此,對於可重複的例子。

library(raster) 
library(sp) 

提出三點rasters並結合成一個stack

set.seed(123) 
r1 <- raster(nrows=10, ncols=10) 
r1 <- setValues(r1, sample(c(1:50), 100, replace = T)) 

r2 <- raster(nrows=10, ncols=10) 
r2 <- setValues(r2, sample(c(40:50), 100, replace = T)) 

r3 <- raster(nrows=10, ncols=10) 
r3 <- setValues(r3, sample(c(50:55), 100, replace = T)) 

Stack <- stack(r1, r2, r3) 
nlayers(Stack) 

充分利用SpatialPoints

Pts <- SpatialPoints(data.frame(x = sample(-175:175, 50, replace = T), 
           y = sample(-75:75, 50, replace = T))) 

情節的柵格之一,點

plot(Stack[[2]]) 
plot(Pts, add = T) 

enter image description here

現在,是否可以更改每個包含至少一個指向NA的單元格的值?在stack上執行此操作會很好。

回答

4

我會用extract(),它可以被要求返回柵格單元的單元數目超過每個點在於:

​​

enter image description here

另外,可以使用rasterize(),雖然它對於非常大的柵格,可能(?)會變慢:

ii <- !is.na(rasterize(Pts, Stack)) 
Stack[ii] <- NA 
+0

很酷。很有幫助。是否有一個原因是你在'Stack'中編入了層'[[1]]'的索引?測試一下,似乎你可以一次對整個堆棧執行操作。例如'ii < - unique(extract(Stack,Pts,cellnumbers = TRUE)[,「cells」])' 'Stack [ii] < - NA' –

+0

@ B.Davis沒有很好的理由,事實上, unique()'也不需要。我已編輯帖子以相應地簡化。 –

相關問題