r - 在R中按高程多边形过滤空间数据帧

Question

我有一个阿尔卑斯山的数字高程模型（包含 x、y、z）。我想提取代表两个海拔高度（<=1000，>1000）的多边形。需要将多边形应用于具有不同 x、y 点 (x,y,snowdepth) 的更精细的阿尔卑斯山配准。必须对 1000m 以上的所有点进行第二次注册过滤。所以。如果我可以存储多边形并将其用于不同的空间数据帧，那就太好了。我的问题是：如何捕获表示 1000m 以上高程的多边形并过滤 R 中的空间数据框。

#packages
packages <- c("RCurl", "RColorBrewer","ggmap","rgeos","fields","dismo","rgdal", "deldir", "dplyr","tidyr","ggplot2","contoureR",
              "maptools","raster","gstat", "magick","lubridate",  "ggplot2","viridis","scales","inlmisc")

has_available   <- packages %in% rownames(installed.packages())
if(any(!has_available)) install.packages(packages[!has_available])

lapply(packages,library,character.only = TRUE)

for(pkg in packages) {
  library(pkg, character.only = TRUE)
}

geo_proj = "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
alps <- shapefile("Alpine_Convention_Perimeter_2018.shp") 

proj4string(alps)
bbalps <- bbox(alps)
alps <- spTransform(alps,geo_proj)

r <- raster(alps,ncols=990,nrows=990) 

srtm <- getData('SRTM', lon=15, lat=47)
srtm2 <- getData('SRTM', lon=12, lat=47)
srtm3 <- getData('SRTM', lon=9, lat=47)
srtm4 <- getData('SRTM', lon=6, lat=45)

#Mosaic/merge srtm tiles
srtmmosaic <- mosaic(srtm, srtm2, srtm3, srtm4, fun=mean)
srtmmosaic <- spTransform(srtmmosaic,geo_proj)
rm(srtm, srtm2, srtm3, srtm4)

rst0 <- projectRaster(srtmmosaic, r, crs=geo_proj)

#set the z-factor, which increases contrast
mosaic <- rst0 * 10 

#then create the hillshade
slope <- terrain(mosaic, opt="slope", unit='radians')
aspect <- terrain(mosaic, opt="aspect", unit='radians')

hillshade <- hillShade(slope, aspect, angle=45, direction=315)
hillshade <- crop(hillshade, extent(alps)) 
hillshade <- mask(hillshade, alps) 

slope <- crop(slope, extent(alps)) 
slope <- mask(slope, alps) 

rst0 <- crop(rst0, extent(alps)) 
rst0 <- mask(rst0,alps) 

#elevation
elev.df <- rasterToPoints(rst0)
elev.df <-  data.frame(elev.df)
colnames(elev.df) <- c("lon", "lat", "elevation")

#ggplot makes me turn the raster into points
hills.df <- rasterToPoints(hillshade)
hills.df <-  data.frame(hills.df)
colnames(hills.df) <- c("lon", "lat", "hills")

#merging the slope shade with the hillshade
slope.df <- rasterToPoints(slope)
slope.df <-  data.frame(slope.df)
colnames(slope.df) <- c("lon", "lat", "slope")
slope.df$slope <- 1- slope.df$slope #invert the scale so that more slope is darker

#elevation normalised
elev.df$val_norm <- (elev.df[,3]-min(elev.df[,3]))/diff(range(elev.df[,3]))

mnt.df<-hills.df %>%
  left_join(slope.df) %>%
  full_join(elev.df)

df = getContourLines(mnt.df, binwidth=1000)
    
ggplot(df,aes(x,y,group=Group,colour=z)) + geom_path()

标签： rgisgeospatial