1

我正在开发一个绘图函数rasterVis::levelplot,用户可以根据它传递一个光栅对象,或者一个光栅对象和一个sf多边形对象。

该函数相当复杂,但显示问题的最小子集如下:

library(sf)
library(raster)
library(rasterVis)

myplot <- function(in_rast, in_poly = NULL) {
  rastplot <- rasterVis::levelplot(in_rast, margin = FALSE)
  polyplot <- layer(sp::sp.polygons(in_poly))
  print(rastplot + polyplot)
}

问题是我在测试时看到了一些奇怪的(对我来说)结果。让我们定义一些虚拟数据 - 一个 1000x1000 栅格和一个sf带有四个多边形的 POYGON 对象,这些多边形将栅格分开 - :

in_rast  <- raster(matrix(nrow = 1000, ncol = 1000))
in_rast  <- setValues(in_rast, seq(1:1000000))

my_poly  <- structure(list(cell_id = 1:4, geometry = structure(list(structure(list(
    structure(c(0, 0.5, 0.5, 0, 0, 0, 0, 0.5, 0.5, 0), .Dim = c(5L, 
    2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
    structure(c(0.5, 1, 1, 0.5, 0.5, 0, 0, 0.5, 0.5, 0), .Dim = c(5L, 
    2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
    structure(c(0, 0.5, 0.5, 0, 0, 0.5, 0.5, 1, 1, 0.5), .Dim = c(5L, 
    2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
    structure(c(0.5, 1, 1, 0.5, 0.5, 0.5, 0.5, 1, 1, 0.5), .Dim = c(5L, 
    2L))), class = c("XY", "POLYGON", "sfg"))), n_empty = 0L, class = c("sfc_POLYGON", 
"sfc"), precision = 0, crs = structure(list(epsg = NA_integer_, 
    proj4string = NA_character_), .Names = c("epsg", "proj4string"
), class = "crs"), bbox = structure(c(0, 0, 1, 1), .Names = c("xmin", 
"ymin", "xmax", "ymax")))), .Names = c("cell_id", "geometry"), row.names = c(NA, 
4L), class = c("sf", "data.frame"), sf_column = "geometry", 
agr = structure(NA_integer_, class = "factor", .Label = c("constant", 
"aggregate", "identity"), .Names = "cell_id"))

并测试功能。从理论上讲,我认为这应该可行:

my_poly  <- as(my_poly, "Spatial") # convert to spatial
myplot(in_rast, in_poly = my_poly)

但我得到:

在此处输入图像描述

这样做:

in_poly <- my_poly
in_poly <- as(in_poly, "Spatial")
myplot(in_rast, in_poly = in_poly)

仍然失败,但结果不同:

在此处输入图像描述

我发现让它工作的唯一方法是从一开始就为多边形对象赋予我在函数中使用的相同名称(即,in_poly):

in_poly <- structure(list(cell_id = 1:4, geometry = structure(list(structure(list(
    structure(c(0, 0.5, 0.5, 0, 0, 0, 0, 0.5, 0.5, 0), .Dim = c(5L, 
    2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
    structure(c(0.5, 1, 1, 0.5, 0.5, 0, 0, 0.5, 0.5, 0), .Dim = c(5L, 
    2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
    structure(c(0, 0.5, 0.5, 0, 0, 0.5, 0.5, 1, 1, 0.5), .Dim = c(5L, 
    2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
    structure(c(0.5, 1, 1, 0.5, 0.5, 0.5, 0.5, 1, 1, 0.5), .Dim = c(5L, 
    2L))), class = c("XY", "POLYGON", "sfg"))), n_empty = 0L, class = c("sfc_POLYGON", 
"sfc"), precision = 0, crs = structure(list(epsg = NA_integer_, 
    proj4string = NA_character_), .Names = c("epsg", "proj4string"
), class = "crs"), bbox = structure(c(0, 0, 1, 1), .Names = c("xmin", 
"ymin", "xmax", "ymax")))), .Names = c("cell_id", "geometry"), row.names = c(NA, 
4L), class = c("sf", "data.frame"), sf_column = "geometry", 
agr = structure(NA_integer_, class = "factor", .Label = c("constant", 
"aggregate", "identity"), .Names = "cell_id"))


in_poly  <- as(in_poly, "Spatial")
myplot(in_rast, in_poly = in_poly)

在此处输入图像描述

谁能解释这里发生了什么?这显然是(?)一个范围界定问题,但我真的不明白为什么函数会这样!

提前致谢 !

4

1 回答 1

1

的帮助页面latticeExtra::layer解释说:

layer 中使用的评估是非标准的,一开始可能会令人困惑:您通常将变量称为面板函数(x、y 等);您通常可以引用存在于全局环境(工作区)中的对象,但在数据参数中按名称将它们传递给 layer.

layer函数内部使用时,您可以将对象嵌入到列表中并将其传递给data参数:

myplot <- function(in_rast, in_poly = NULL) {
  rastplot <- levelplot(in_rast, margin = FALSE)
  polyplot <- layer(sp.polygons(x), 
                    data = list(x = in_poly))
  print(rastplot + polyplot)
}

现在该函数产生所需的结果:

myplot(in_rast, in_poly = my_poly)

在此处输入图像描述

于 2017-07-21T06:41:46.413 回答