# 加载运行环境
library(tidyverse)
library(sf)
library(raster)
library(rgdal)
library(glue)
library(data.table)
library(ggsn)
library(cowplot)
## read dem file
myfile <- "00DataSet/dem_250m/dem_250m/w001001.adf"
dem <- new("GDALReadOnlyDataset", myfile) %>% asSGDF_GROD
crs.geo <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84")
## get the boundary of Three River Plain
sjpy <-
read_sf("00DataSet/Province/CN-sheng-A.shp",
options = "ENCODING=gb2312",
stringsAsFactors = FALSE) %>%
filter(name == "辽宁") %>%
st_transform(crs = crs.geo)
bdf <-sjpy %>% getElement("geometry") %>%
st_transform(crs = crs.geo) %>% as_Spatial()
## mask dem using the boundary
dem.sj <- raster(dem)
dem.sj <- dem.sj %>% crop(bdf) %>%
raster::aggregate(fact = 0.01 / res(d
利用R绘制指定区域的DEM图
最新推荐文章于 2024-02-17 21:21:39 发布
![](https://img-home.csdnimg.cn/images/20240711042549.png)