通过 R 中的高程多边形过滤空间数据框

编程问答 2022-05-10

问题描述

我有一个阿尔卑斯山的数字高程模型（包含 x、y、z）。我想提取代表两个高程级别（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","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"
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,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","hills")

#merging the slope shade with the hillshade
slope.df <- rasterToPoints(slope)
slope.df <-  data.frame(slope.df)
colnames(slope.df) <- c("lon","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,group=Group,colour=z)) + geom_path()

解决方法

暂无找到可以解决该程序问题的有效方法，小编努力寻找整理中！

如果你已经找到好的解决方法，欢迎将解决方案带上本链接一起发送给小编。

小编邮箱:dio#foxmail.com (将#修改为@）

geospatialgisrr