使用 R 模拟海平面变化

Question

Dropbox linked files 我正在尝试使用 R 模拟海平面变化，但不知道该怎么做。类似的期望结果。 Similar desired result 我看过一篇文章说有人使用了 contour() 但不知道如何使用这个函数。我还看到很多人使用 ggplot2，但我不确定如何使用我的数据集来实现数据框。

我有一个 dtm，它是 5 个 asc tiles 的组合版本。使用下面的代码。

setwd()
f <-list.files(pattern = ".asc")  
r <- lapply(f, raster) 
x <- do.call("merge",r) 
writeRaster(x,"DTM_combine.asc", overwrite=TRUE)  
library(rgdal)
r = raster("DTM_combine.asc")
plot(r)

然后我使用这个对栅格进行了重新分类。

image(r,zlim=c(0,70), main="DEM Findhornbay", col=col) m=c(0,5,1,5,10,2,10,15,3,15,20,4,20,25,5,25,30,6,30,35,7,35,40,8,40,45,9,45,50,10,50,55,11,55,60,12,60,65,13,65,70,14)
mat=matrix(m,ncol=3,byrow=TRUE)
r=reclassify(r,mat)
rcat=reclassify(r,mat)
plot(rcat)
col <-terrain.colors(14)
plot(rcat)
brk <-c(1,2,3,4,5,6,7,8,9,10,11,12,13,14)
plot(rcat, col=col, breaks=brk, main="Recalssed DEM Findhorn Bay")
plot(r, col=col, breaks=brk, main="Recalssed DEM Findhorn Bay")

产生了这个。我不确定这是否是正确的工作方式，但下一个 objective 是模拟海平面变化。

plot3D(r)

Answer 1

栅格包和底图

要使用光栅和光栅包提供的功能绘制等高线，您可以按如下方式进行：

首先加载包，声明您的 crs，加载您的栅格并分配您的 crs：

library(raster)
library(magrittr)

#We define the crs of the project
crs <- "+proj=utm +datum=WGS84 +no_defs"

#read your raster DTM_combine.asc and assign it the projection
r <- raster::raster("PathToFolder/DTM_combine.asc")
raster::projection(r) <- crs

就我们使用数字高程模型而言，我假设高程单位为米。我们可以将海拔矢量 niveaux 定义为每 5 米介于 0 和最大海拔值之间的序列：

niveaux <- (raster::extract(r, raster::extent(r)) %>% 
              range(finite = TRUE))[2] %>% 
                seq(0,.,by=5)

然后你可以绘制光栅：

raster::plot(r)

然后添加轮廓。假设海平面上升 5 米，这里我添加了蓝色的海岸线（海拔 =0）和后退的海岸线（红色）。

raster::contour(r, add=TRUE, col="blue", lwd=0.2, levels=niveaux[1], drawlabels=FALSE)
raster::contour(r, add=TRUE, col="red", lwd=0.2, levels=niveaux[2], drawlabels=FALSE)

或者，如果按如下方式绘制等高线，则可以只使用等高线：

raster::contour(r, col="blue", lwd=0.2, levels=niveaux[1], drawlabels=FALSE)
raster::contour(r, add=TRUE, col="red", lwd=0.2, levels=niveaux[2], drawlabels=FALSE)

ggplot2 方式

首先，将栅格转换为数据框：

r.df <- raster::as.data.frame(r, xy = TRUE)

然后用 geom_raster 绘制光栅，用 geom_contour 绘制等高线，方法与之前相同：

  ggplot(data = r.df, aes(x = x, y = y)) +
  geom_raster(aes(fill=DTM_combine)) +
  geom_contour(aes(z=DTM_combine, color = "Current"), breaks = niveaux[1], size=0.5) + 
  geom_contour(aes(z=DTM_combine, color = "New"), breaks = niveaux[2], size=0.5) +
  scale_fill_gradientn(colours = terrain.colors(10, rev=TRUE), na.value = "transparent") +
  theme_void() +
  labs(fill="Elevation", color="Coastal lines")

生成的地图：

如果您想避开左上角的白色方块，您可以将 scale_fill_gradientn() 中的 na.value 修改为 #f2f2f2 例如。

*奖励：如果你想填充海岸线之间的 space 你可以看看 [=28] 中的函数 geom_contour_fill() =]包。

希望对您有所帮助

使用 R 模拟海平面变化

Simulate sea level change using R

r

simulate

ggplot2