将 XY 点添加到 levelplot 生成的栅格地图答案

【问题标题】：Add XY points to raster map generated by levelplot将 XY 点添加到 levelplot 生成的栅格地图
【发布时间】：2015-02-19 01:21:27
【问题描述】：

我有使用 R 中的 raster 包生成的 raster 地图。这些栅格图层可以使用 rasterVis 包的 levelplot 函数进行可视化：

levelplot(rasterstack, layout=c(1, 2), 
          col.regions=colorRampPalette(c('darkred', 'red3', 'orange2', 'orange', 
                                         'yellow', 'lightskyblue', 'steelblue3', 
                                         'royalblue3', 'darkblue')))

现在，我想将一些由 xy 坐标定义的 z 值添加到 levelplot 地图。包含 z 值的数据框有 4 列。第 1 列和第 2 列包含 x 和 y 坐标，第 3 列包含布局 (1, 1) 中地图 1 的 z 值，第 4 列包含布局(1, 2) 的 z 值。

应该添加每个地图的点，如果 z pch=2，如果 z > 0.05，pch=3。

我在网上搜索并找到了 Ripley 的解决方案，但在我的情况下它不起作用：

levelplot(rcp852, xlab = "", ylab = "",
          panel = function(x, y, subscripts, ...) {
            panel.levelplot(x, y, subscripts, ...)
            panel.xyplot(topo$x,topo$y, cex = 0.5, col = 1)
          }
)

我尝试了许多其他选项，但这些点与通过levelplot 生成的地图不对齐。

【问题讨论】：

标签： r points levelplot

【解决方案1】：

layer 这个很方便：

library(raster)
library(rasterVis)
library(sp)

s <- stack(replicate(2, raster(matrix(runif(100), 10))))
xy <- data.frame(coordinates(sampleRandom(s, 10, sp=TRUE)),
                z1=runif(10), z2=runif(10))
coordinates(xy) <- ~x+y

levelplot(s, margin=FALSE, at=seq(0, 1, 0.05)) + 
  layer(sp.points(xy, pch=ifelse(xy$z1 < 0.5, 2, 3), cex=2, col=1), columns=1) +
  layer(sp.points(xy, pch=ifelse(xy$z2 < 0.5, 2, 3), cex=2, col=1), columns=2)

请注意，layer 的 columns 参数（rows 也存在）指定您要将图层添加到哪个面板。

【讨论】：

完美！你刚刚为我节省了世界时间。
它现在应该可以再次使用@JerryN。感谢您的提醒。

【解决方案2】：

所以我有一个带有 xy 坐标和许多 z 列的数据框。感谢@jbaums，这是将积分添加到我的地图的最终答案：

s <- stack(raster1,raster2)
coordinates(SITES...TTEST) <- ~x+y # SpatialPointsDataFrame
levelplot(s, layout=c(1, 2), 
          col.regions=colorRampPalette(c("darkred", "red3", 
                                         "orange", "yellow", "lightskyblue", "royalblue3", 
                                         "darkblue")),
          at=seq(floor(39.15945) ,ceiling(51.85068), length.out=30), 
          par.strip.text=list(cex=0),scales=list(alternating=FALSE))+
  layer(sp.points(SITES...TTEST, pch=ifelse(SITES...TTEST$Precip_DJF6 < 0.05, 2, 3), cex=2, col=1), rows=1) +
  layer(sp.points(SITES...TTEST, pch=ifelse(SITES...TTEST$Precip_DJF6 < 0.05, 2, 3), cex=2, col=1), rows=2)

【讨论】：