1## Convert 'SpatialLines*' object to spatstat 'linnet' object
2##
3## For 'SpatialLinesDataFrame', the data columns are copied
4## to the network as marks associated with the network segments.
5##
6## If fuse=TRUE, the code searches for pairs of points with the same (x,y)
7##               coordinates that occur in different polylines,
8##               and merges them together as identical vertices of the network.
9##
10##    Last edit: 2020/04/20 Adrian Baddeley
11
12if (!isClass("linnet"))
13	setClass("linnet")
14
15as.linnet.SpatialLines <- function(X, ..., fuse=TRUE) {
16  if (!is.na(sp::is.projected(X)) && !sp::is.projected(X))
17    stop("Only projected coordinates may be converted to spatstat class objects")
18  check_spatstat("spatstat.geom")
19  check_spatstat("spatstat.linnet")
20  # extract bounding box to use as window
21  bb <- bbox(X)
22  BB <- spatstat.geom::owin(bb[1,], bb[2,])
23  #
24  n <- length(X)
25  xx <- yy <- numeric(0)
26  ii <- jj <- integer(0)
27  if(n > 0) {
28    # coordinates of all vertices
29    crdlists <- coordinates(X)
30    rowcounts <- unname(lapply(crdlists, function(x) sapply(x, nrow)))
31    colcounts <- unname(lapply(crdlists, function(x) sapply(x, ncol)))
32    if(any(unlist(colcounts) != 2)) stop("Coordinates should be 2-column matrices", call.=FALSE)
33    # 'rbind' all the matrices of coordinates
34    xy <- unlist(lapply(crdlists, function(x) lapply(x, t)))
35    xy <- matrix(xy, ncol=2, byrow=TRUE)
36    xx <- xy[,1]
37    yy <- xy[,2]
38    # construct indices for each level of list in original data
39    ii <- rep(seq_along(crdlists), sapply(rowcounts, sum))
40    jj <- unlist(lapply(rowcounts, function(x) rep(seq_along(x), as.integer(x))))
41    # check for *repeated* vertices within the same line
42    rpt <- c(FALSE, (diff(xx) == 0) & (diff(yy) == 0) & (diff(ii) == 0) & (diff(jj) == 0))
43    if(any(rpt)) {
44      warning("Repeated vertices (on the same line) were removed", call.=FALSE)
45      retain <- !rpt
46      xx <- xx[retain]
47      yy <- yy[retain]
48      ii <- ii[retain]
49      jj <- jj[retain]
50    }
51  }
52  # extract vertices
53  V <- spatstat.geom::ppp(xx, yy, window=BB, check=!fuse)
54  nV <- length(xx)
55  # join them
56  edges <- NULL
57  iii <- jjj <- integer(0)
58  if(nV > 1) {
59    seqn <- seq_len(nV)
60    from <- seqn[-nV]
61    to   <- seqn[-1]
62    ok   <- diff(ii) == 0 & diff(jj) == 0
63    from <- from[ok]
64    to   <- to[ok]
65    iii  <- ii[c(ok, FALSE)] # indices backward
66    jjj  <- jj[c(ok, FALSE)]
67    if(fuse) {
68      umap <- spatstat.geom::uniquemap(V)
69      retain <- (umap == seq_along(umap))
70      V <- V[retain]
71      renumber <- cumsum(retain)
72      from <- renumber[umap[from]]
73      to   <- renumber[umap[to]]
74    }
75    edges <- cbind(from, to)
76  }
77  if(!is.null(edges)) {
78    up <- (from < to)
79    ee <- cbind(ifelse(up, from , to), ifelse(up, to, from))
80    if(anyDuplicated(ee)) {
81      u <- !duplicated(ee)
82      from <- from[u]
83      to   <- to[u]
84      iii  <- iii[u]
85      jjj  <- jjj[u]
86    }
87  }
88  result <- spatstat.linnet::linnet(vertices=V, edges = edges, sparse=TRUE)
89  if(spatstat.geom::nsegments(result) == length(iii)) {
90    df <- data.frame(LinesIndex=iii, LineIndex=jjj)
91    if(.hasSlot(X, "data")) {
92      DF <- slot(X, "data")
93      df <- cbind(DF[iii,,drop=FALSE], df)
94    }
95    spatstat.geom::marks(result$lines) <- df
96  } else warning("Internal error: could not map data frame to lines")
97  return(result)
98}
99
100setAs("SpatialLines", "linnet", function(from) as.linnet.SpatialLines(from))
101
102setAs("SpatialLinesDataFrame", "linnet", function(from) as.linnet.SpatialLines(from))
103
104