【问题标题】:Smoothing map with elevation levels; R with tmap用高程平滑地图; R 与 tmap
【发布时间】:2020-06-29 21:48:30
【问题描述】:

我正在为一份科学手稿制作地图。该地图应显示墨西哥的基本地图,带有采样点并显示海拔高度颜色。为此,我将 R 与包 tmap 一起使用。这是我的代码:

library(sf)
library(raster)
library(dplyr)
library(spData)
library(spDataLarge)
library(tmap)    # for static and interactive maps

data(World)
data("metro")
data("land")

gpsdata <- read.csv('gps.csv')
gpsdata <- as.matrix(gpsdata)
gpspoints <- SpatialPoints(gpsdata)
plot(gpspoints)

gpsnames <- read.csv('gps_names.csv')

spdf = SpatialPointsDataFrame(gpsdata,gpsnames)

proj4string(gpspoints) <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
proj4string(spdf) <- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")

pdf('fig-sampling.pdf')
tm_shape(World[World$name=="Mexico", ]) + 
  tm_polygons() +
tm_shape(land) +
  tm_raster("elevation", breaks=c(-Inf, 250, 500, 1000, 1500, 2000, 2500, 3000, 4000, Inf),  
            palette = terrain.colors(9), title="Elevation") +
tm_shape(World) +
  tm_borders("black", lwd = 1) +
  tm_text("name") +
tm_layout(legend.position = c("right","center"),
          legend.bg.color = "lightblue", bg.color = 'lightblue') +
tm_shape(metro) +
  tm_bubbles("pop2020", title.size = "Population") +
  tm_text("name", size = "pop2010", legend.size.show = FALSE, root=8, size.lowerbound = 0.7, just = 'right', xmod = -0.5, ymod = 0)+
  tm_grid(projection="longlat", labels.size = .5) +
tm_shape(spdf) +
  tm_bubbles(col = 'blue', size = 0.3) +
  tm_text("Location", size = 0.8, just = 'left', xmod = 0.3)
dev.off()

不幸的是,高程图层以非常大的像素绘制:

有什么方法可以平滑(我已经尝试过 style = 'cont' 用于高程栅格图层)?

【问题讨论】:

    标签: r maps scientific-computing geography tmap


    【解决方案1】:

    最简单的方法是使用 raster 包中的 disaggregate 对其进行插值。

    这会使栅格放大给定因子,并插入影响平滑度的像素值。但是对于您的计算机(它是我的)来说,整个世界的更大栅格可能太多了,所以首先裁剪到墨西哥。

    世界数据与土地位于不同的坐标系中,因此您必须将其转换为作物的土地坐标:

    > land = crop(land, st_transform(World[World$name=="Mexico",],crs(land)))
    > plot(land)
    

    这应该只显示墨西哥,但仍然是块状的。

    > land = disaggregate(land, 4, "bilinear")
    > plot(land)
    

    这应该显示墨西哥,但平滑。请注意,插值似乎有伪像,不应该用于分析,它们只是看起来更漂亮一些。

    如果您想进行分析,可以通过 SRTM 数据集和其他全球 DEM 数据集获得 30 米分辨率的高程数据。

    【讨论】:

    • 非常感谢@Spacedman!我稍微修改了您的代码以包含邻近的土地区域并增加了“平滑”的参数:landsel &lt;- extent(-120,-80,12,34)land = crop(land, landsel)land = disaggregate(land, 10, "bilinear")plot(land$elevation, col = terrain.colors(9)) 不完美(政治边界和地形图有一些差异),但对于字段/采样图。
    • 您可以尝试使用raster::mask(land, World[....]) 之类的方法将平滑栅格遮盖到World.... 的陆地/海洋边界(需要st_project World 同上)
    猜你喜欢
    • 1970-01-01
    • 2021-11-03
    • 2018-06-21
    • 2015-05-07
    • 2020-04-06
    • 2021-06-27
    • 1970-01-01
    • 2013-06-06
    • 2017-06-01
    相关资源
    最近更新 更多