如果在其他栅格中满足条件,则尝试使用哪个函数从栅格中提取数据
trying to use which function to pull data from a raster if condition satisfied in other raster
我有两个大小相同且包含来自同一位置但数据类型不同的数据的栅格(一个栅格具有坡度数据,另一个栅格具有坡向数据)。我希望能够一次查看一个方面的坡度数据,所以我试图创建一个设置(也许是 if/else 语句?)我说“如果(方面条件)在一个栅格中得到满足,坡度数据将从另一个栅格中的同一像素中提取。
#I have a slope and an aspect raster that i pulled
library(raster)
library(rgdal)
library(sp)
aspect <- raster("geotiff name here")
slope <- raster("geotiff name here")
#Looking at the north aspect (between 0-22.5 degrees or 337.5-360 degrees)
#First I am setting the pixels in the aspect raster that correspond to north
#equal to 1, and the values that don't = 0
aspect[aspect >= 0 & aspect <= 22.5] <- 1
aspect[aspect >= 337.5 & aspect <= 360] <- 1
aspect[aspect > 22.5 & aspect < 337.5] <- 0
#Here i am saving the indices of the raster that face north to a new one
north <- which(aspect == 1, cells = true)
然后我只想从坡度栅格的像素中读取数据,这些像素从坡向栅格中分配了 TRUE 值,但这就是我被难住的地方!我最近开始使用 R,所以可能有一种简单的方法可以做到这一点,但我很想念它,我们将不胜感激。非常感谢!
您不需要将 1 转换为 TRUE,因为 R 会自动执行此操作。试试这个代码:
#create a data frame
data <- data.frame(aspect=aspect, slope=slope)
#create a 'north' column and populate with 1
data$north <- 1
#those that don't meet the north criteria are converted to 0
data$north[data$aspect > 22.5 & data$aspect < 337.5] <- 0
#report the 'slope' values where north=1
data$slope[data$north == 1]
始终包含示例数据(请参阅帮助文件以获取灵感,此处来自 ?raster::terrain
)
library(raster)
x <- getData('alt', country='CHE')
aspect <- terrain(x, 'aspect', unit='degrees')
slope <- terrain(x, 'slope', unit='degrees')
这是一种更好的重分类方法:
m <- matrix(c(0,22.5,1,22.5,337.50,0,337.5,360,1), ncol=3, byrow=TRUE)
aspectcls <- reclassify(aspect, m)
获取坡度数据,其中 aspectcls != 0
nslope <- mask(slope, aspectcls, maskvalue=0)
获取值
v <- values(nslope)
boxplot(v)
你也可以
crosstab(aspectcls, slope)
我不推荐你走的路,但如果你走的话,你可以
cells <- Which(aspectcls, cells=T)
vv <- slope[cells]
boxplot(vv)
我有两个大小相同且包含来自同一位置但数据类型不同的数据的栅格(一个栅格具有坡度数据,另一个栅格具有坡向数据)。我希望能够一次查看一个方面的坡度数据,所以我试图创建一个设置(也许是 if/else 语句?)我说“如果(方面条件)在一个栅格中得到满足,坡度数据将从另一个栅格中的同一像素中提取。
#I have a slope and an aspect raster that i pulled
library(raster)
library(rgdal)
library(sp)
aspect <- raster("geotiff name here")
slope <- raster("geotiff name here")
#Looking at the north aspect (between 0-22.5 degrees or 337.5-360 degrees)
#First I am setting the pixels in the aspect raster that correspond to north
#equal to 1, and the values that don't = 0
aspect[aspect >= 0 & aspect <= 22.5] <- 1
aspect[aspect >= 337.5 & aspect <= 360] <- 1
aspect[aspect > 22.5 & aspect < 337.5] <- 0
#Here i am saving the indices of the raster that face north to a new one
north <- which(aspect == 1, cells = true)
然后我只想从坡度栅格的像素中读取数据,这些像素从坡向栅格中分配了 TRUE 值,但这就是我被难住的地方!我最近开始使用 R,所以可能有一种简单的方法可以做到这一点,但我很想念它,我们将不胜感激。非常感谢!
您不需要将 1 转换为 TRUE,因为 R 会自动执行此操作。试试这个代码:
#create a data frame
data <- data.frame(aspect=aspect, slope=slope)
#create a 'north' column and populate with 1
data$north <- 1
#those that don't meet the north criteria are converted to 0
data$north[data$aspect > 22.5 & data$aspect < 337.5] <- 0
#report the 'slope' values where north=1
data$slope[data$north == 1]
始终包含示例数据(请参阅帮助文件以获取灵感,此处来自 ?raster::terrain
)
library(raster)
x <- getData('alt', country='CHE')
aspect <- terrain(x, 'aspect', unit='degrees')
slope <- terrain(x, 'slope', unit='degrees')
这是一种更好的重分类方法:
m <- matrix(c(0,22.5,1,22.5,337.50,0,337.5,360,1), ncol=3, byrow=TRUE)
aspectcls <- reclassify(aspect, m)
获取坡度数据,其中 aspectcls != 0
nslope <- mask(slope, aspectcls, maskvalue=0)
获取值
v <- values(nslope)
boxplot(v)
你也可以
crosstab(aspectcls, slope)
我不推荐你走的路,但如果你走的话,你可以
cells <- Which(aspectcls, cells=T)
vv <- slope[cells]
boxplot(vv)