如何检查栅格和多边形是否重叠?
How to check if raster and polygon are overlapping?
我有一个栅格,我想检查具有值的像元是否与我感兴趣的区域重叠,这是一个空间多边形文件。
有快速的方法吗?
我的光栅看起来像这样:
我尝试了什么:
r <- raster("C:/user/test.jp2")
studyarea <- readOGR("C:/user/boundary.shp")
dummi <- as(extent(r), "SpatialPolygons")
if (gContainsProperly(studyarea, dummi)) {
print("Raster extent is fully within the Area of interest")
}
else if (gIntersects(studyarea, dummi)) {
print("Raster extent is not fully within the Area of interest")
} else {
print("Raster extent is fully outside the Area of interest")
}
as(extent(r), "SpatialPolygons")
在值和 NA 中创建一个多边形。但我只想检查值部分和研究区域是否重叠
栅格和研究区域的数据示例:
matrix <- matrix(c(NA,NA,NA,1,NA,NA,2,1,NA,1,3,4,2,3,4,1), nrow=4)
r <- raster(matrix)
extent(r) <- c(399960, 509760, 5290200, 5400000)
crs(r) <- "+proj=utm +zone=32 +datum=WGS84 +units=m +no_defs"
writeRaster(r, paste(getwd(),"filename.jp2"))
matrix <- matrix(c(1,1,NA,NA,1,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA), nrow=4)
rpoly<- raster(matrix)
extent(rpoly) <- c(399960, 509760, 5290200, 5400000)
crs(rpoly) <- "+proj=utm +zone=32 +datum=WGS84 +units=m +no_defs"
studyarea <- rasterToPolygons(rpoly)
谢谢
所以这是使用数据示例的方法。
想法是将具有值的栅格像元转换为单个多边形,并将该多边形与您的研究区域多边形进行比较。
sf
提供了多种谓词来比较几何。
matrix <- matrix(c(NA,NA,NA,1,NA,NA,2,1,NA,1,3,4,2,3,4,1), nrow=4)
r <- raster::raster(matrix)
raster::extent(r) <- c(399960, 509760, 5290200, 5400000)
raster::crs(r) <- "+proj=utm +zone=32 +datum=WGS84 +units=m +no_defs"
r <- sf::st_as_sf(raster::rasterToPolygons(r)) # This will create a sfc polygon dataframe from your cells with values
r <- sf::st_union(r) # This will combine all cells with values to a single polygon
matrix <- matrix(c(1,1,NA,NA,1,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA), nrow=4)
rpoly<- raster::raster(matrix)
raster::extent(rpoly) <- c(399960, 509760, 5290200, 5400000)
raster::crs(rpoly) <- "+proj=utm +zone=32 +datum=WGS84 +units=m +no_defs"
studyarea <- sf::st_as_sf(raster::rasterToPolygons(rpoly)) # This will create a sfc polygon dataframe from your cells with values
studyarea <- sf::st_union(studyarea) # This will combine all cells with values to a single polygon
sf::st_contains(r, studyarea, sparse = FALSE) # This checks if studyarea is within r, see ?sf::st_contains for other predicates you might want (touches, crosses, etc.)
我找到了一个更快的解决方案。可以通过多边形裁剪栅格,然后通过多边形屏蔽栅格。如果多边形区域中没有栅格值。栅格将具有最小值 = 0 和最大值 = 0。
创建数据
matrix <- matrix(c(NA,NA,NA,1,NA,NA,2,1,NA,1,3,4,2,3,4,1), nrow=4)
r <- raster(matrix)
extent(r) <- c(399960, 509760, 5290200, 5400000)
crs(r) <- "+proj=utm +zone=32 +datum=WGS84 +units=m +no_defs"
writeRaster(r, paste(getwd(),"filename.jp2"))
matrix <- matrix(c(1,1,NA,NA,1,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA), nrow=4)
rpoly<- raster(matrix)
extent(rpoly) <- c(399960, 509760, 5290200, 5400000)
crs(rpoly) <- "+proj=utm +zone=32 +datum=WGS84 +units=m +no_defs"
studyarea <- rasterToPolygons(rpoly)
解决方案
# Reproject. Make sure both have the same projection
studyarea <- spTransform(studyarea, projection(r))
# Crop and clip raster with study area polygon
r_crop <- crop(r, studyarea)
r_mask <- mask(r_crop, studyarea, updatevalue = NA, updateNA= TRUE)
if (r_mask@data@min == 0 && r_mask@data@max == 0){
print("raster and studyarea do not overlap))
}
else {
print("Raster and studyarea do overlap")
}
我有一个栅格,我想检查具有值的像元是否与我感兴趣的区域重叠,这是一个空间多边形文件。
有快速的方法吗?
我的光栅看起来像这样:
我尝试了什么:
r <- raster("C:/user/test.jp2")
studyarea <- readOGR("C:/user/boundary.shp")
dummi <- as(extent(r), "SpatialPolygons")
if (gContainsProperly(studyarea, dummi)) {
print("Raster extent is fully within the Area of interest")
}
else if (gIntersects(studyarea, dummi)) {
print("Raster extent is not fully within the Area of interest")
} else {
print("Raster extent is fully outside the Area of interest")
}
as(extent(r), "SpatialPolygons")
在值和 NA 中创建一个多边形。但我只想检查值部分和研究区域是否重叠
栅格和研究区域的数据示例:
matrix <- matrix(c(NA,NA,NA,1,NA,NA,2,1,NA,1,3,4,2,3,4,1), nrow=4)
r <- raster(matrix)
extent(r) <- c(399960, 509760, 5290200, 5400000)
crs(r) <- "+proj=utm +zone=32 +datum=WGS84 +units=m +no_defs"
writeRaster(r, paste(getwd(),"filename.jp2"))
matrix <- matrix(c(1,1,NA,NA,1,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA), nrow=4)
rpoly<- raster(matrix)
extent(rpoly) <- c(399960, 509760, 5290200, 5400000)
crs(rpoly) <- "+proj=utm +zone=32 +datum=WGS84 +units=m +no_defs"
studyarea <- rasterToPolygons(rpoly)
谢谢
所以这是使用数据示例的方法。
想法是将具有值的栅格像元转换为单个多边形,并将该多边形与您的研究区域多边形进行比较。
sf
提供了多种谓词来比较几何。
matrix <- matrix(c(NA,NA,NA,1,NA,NA,2,1,NA,1,3,4,2,3,4,1), nrow=4)
r <- raster::raster(matrix)
raster::extent(r) <- c(399960, 509760, 5290200, 5400000)
raster::crs(r) <- "+proj=utm +zone=32 +datum=WGS84 +units=m +no_defs"
r <- sf::st_as_sf(raster::rasterToPolygons(r)) # This will create a sfc polygon dataframe from your cells with values
r <- sf::st_union(r) # This will combine all cells with values to a single polygon
matrix <- matrix(c(1,1,NA,NA,1,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA), nrow=4)
rpoly<- raster::raster(matrix)
raster::extent(rpoly) <- c(399960, 509760, 5290200, 5400000)
raster::crs(rpoly) <- "+proj=utm +zone=32 +datum=WGS84 +units=m +no_defs"
studyarea <- sf::st_as_sf(raster::rasterToPolygons(rpoly)) # This will create a sfc polygon dataframe from your cells with values
studyarea <- sf::st_union(studyarea) # This will combine all cells with values to a single polygon
sf::st_contains(r, studyarea, sparse = FALSE) # This checks if studyarea is within r, see ?sf::st_contains for other predicates you might want (touches, crosses, etc.)
我找到了一个更快的解决方案。可以通过多边形裁剪栅格,然后通过多边形屏蔽栅格。如果多边形区域中没有栅格值。栅格将具有最小值 = 0 和最大值 = 0。
创建数据
matrix <- matrix(c(NA,NA,NA,1,NA,NA,2,1,NA,1,3,4,2,3,4,1), nrow=4)
r <- raster(matrix)
extent(r) <- c(399960, 509760, 5290200, 5400000)
crs(r) <- "+proj=utm +zone=32 +datum=WGS84 +units=m +no_defs"
writeRaster(r, paste(getwd(),"filename.jp2"))
matrix <- matrix(c(1,1,NA,NA,1,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA), nrow=4)
rpoly<- raster(matrix)
extent(rpoly) <- c(399960, 509760, 5290200, 5400000)
crs(rpoly) <- "+proj=utm +zone=32 +datum=WGS84 +units=m +no_defs"
studyarea <- rasterToPolygons(rpoly)
解决方案
# Reproject. Make sure both have the same projection
studyarea <- spTransform(studyarea, projection(r))
# Crop and clip raster with study area polygon
r_crop <- crop(r, studyarea)
r_mask <- mask(r_crop, studyarea, updatevalue = NA, updateNA= TRUE)
if (r_mask@data@min == 0 && r_mask@data@max == 0){
print("raster and studyarea do not overlap))
}
else {
print("Raster and studyarea do overlap")
}