snapPointsToLines无法保留R中的属性信息

3

最近我发现 snapPointsToLines 存在一个问题,它无法保留空间点数据框的属性。以下是一个示例:

# Generate a spatial line dataframe
l1 = cbind(c(1,2,3),c(3,2,2))
l1a = cbind(l1[,1]+.05,l1[,2]+.05)
l2 = cbind(c(1,2,3),c(1,1.5,1))
Sl1 = Line(l1)
Sl1a = Line(l1a)
Sl2 = Line(l2)
S1 = Lines(list(Sl1, Sl1a), ID="a")
S2 = Lines(list(Sl2), ID="b")
Sl = SpatialLines(list(S1,S2))
df = data.frame(z = c(1,2), row.names=sapply(slot(Sl, "lines"), function(x) slot(x, "ID")))
Sldf = SpatialLinesDataFrame(Sl, data = df)

# Generate a spatial point dataframe
xc = c(1.2,1.5,2.5)
yc = c(1.5,2.2,1.6)
Spoints = SpatialPoints(cbind(xc, yc))
Spdf <- SpatialPointsDataFrame(Spoints, data = data.frame(value = 1:length(Spoints)))

#use the function SpatialPointsDataFrame 
res <- snapPointsToLines(Spdf, Sldf)

res只有“nearest_line_id”和“snap_dist”两个字段,没有我需要的Spdf中的“value”字段。

#use the function SpatialPointsDataFrame with "withAttrs = TRUE" parameter
res <- snapPointsToLines(Spdf, Sldf, withAttrs = TRUE)

报告错误:

"Error in snapPointsToLines(Spdf, Sldf, withAttrs = TRUE) : 
  A SpatialPoints object has no attributes! Please set withAttrs as FALSE."

但是 Spdf 是带属性的空间点数据框。我不知道它出了什么问题。几周前我使用这个函数时,它没有这个问题。

1个回答

2

我认为问题可能是由函数本身引起的。当你查看这个函数的代码时,我们可以看到下面的开头部分。

if (class(points) == "SpatialPoints" && missing(withAttrs)) 
    withAttrs = FALSE
if (class(points) == "SpatialPoints" && withAttrs == TRUE) 
    stop("A SpatialPoints object has no attributes! Please set withAttrs as FALSE.")

有时候,SpatialPointsDataFrame可能被识别为SpatialPoints。因此,该函数将把您的SpatialPointsDataFrame视为SpatialPoints,并不会在函数中保留属性。 您可以按照下面的代码对函数进行一些修改。
snapPointsToLines1 <-  function (points, lines, maxDist = NA, withAttrs = TRUE, idField = NA) 
{
  if (rgeosStatus()) {
    if (!requireNamespace("rgeos", quietly = TRUE)) 
      stop("package rgeos required for snapPointsToLines")
  }
  else stop("rgeos not installed")
  if (is(points, "SpatialPointsDataFrame")==FALSE && missing(withAttrs)) 
    withAttrs = FALSE
  if (is(points, "SpatialPointsDataFrame")==FALSE && withAttrs == TRUE) 
    stop("A SpatialPointsDataFrame object is needed! Please set withAttrs as FALSE.")
  d = rgeos::gDistance(points, lines, byid = TRUE)
  if (!is.na(maxDist)) {
    distToLine <- apply(d, 2, min, na.rm = TRUE)
    validPoints <- distToLine <= maxDist
    distToPoint <- apply(d, 1, min, na.rm = TRUE)
    validLines <- distToPoint <= maxDist
    points <- points[validPoints, ]
    lines = lines[validLines, ]
    d = d[validLines, validPoints, drop = FALSE]
    distToLine <- distToLine[validPoints]
    if (!any(validPoints)) {
      if (is.na(idField)) {
        idCol = character(0)
      }
      else {
        idCol = lines@data[, idField][0]
      }
      newCols = data.frame(nearest_line_id = idCol, snap_dist = numeric(0))
      if (withAttrs) 
        df <- cbind(points@data, newCols)
      else df <- newCols
      res <- SpatialPointsDataFrame(points, data = df, 
                                    proj4string = CRS(proj4string(points)), match.ID = FALSE)
      return(res)
    }
  }
  else {
    distToLine = apply(d, 2, min, na.rm = TRUE)
  }
  nearest_line_index = apply(d, 2, which.min)
  coordsLines = coordinates(lines)
  coordsPoints = coordinates(points)
  mNewCoords = vapply(1:length(points), function(x) nearestPointOnLine(coordsLines[[nearest_line_index[x]]][[1]], 
                                                                       coordsPoints[x, ]), FUN.VALUE = c(0, 0))
  if (!is.na(idField)) {
    nearest_line_id = lines@data[, idField][nearest_line_index]
  }
  else {
    nearest_line_id = sapply(slot(lines, "lines"), 
                             function(i) slot(i, "ID"))[nearest_line_index]
  }
  if (withAttrs) 
    df = cbind(points@data, data.frame(nearest_line_id, snap_dist = distToLine))
  else df = data.frame(nearest_line_id, snap_dist = distToLine, 
                       row.names = names(nearest_line_index))
  SpatialPointsDataFrame(coords = t(mNewCoords), data = df, 
                         proj4string = CRS(proj4string(points)))
}

然后,使用这个新函数snapPointsToLines1,你可以获得你想要的属性。

有人为此提交了错误报告吗?毕竟这是一个相当荒谬的错误 ;)! - Gmichael

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接