daily pastebin goal
6%
SHARE
TWEET

Untitled

a guest Feb 22nd, 2019 63 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. require(rgdal)  
  2.     require(rgeos)
  3.     library(sp)
  4.     library(raster)
  5.     library(spatialEco)
  6.     library(spdep)
  7.  
  8.     Shapefile <- readOGR(".","US_place_2000")
  9.     Shapefile <- sp.na.omit(Shapefile, margin=1) # NAs are present for some reason?
  10.     Touching_List <- gTouches(Shapefile, byid = TRUE, returnDense=FALSE)
  11.     # also tried Touching_List <- poly2nb(Shapefile, queen = FALSE)
  12.  
  13.     all.length.list <- lapply(1:length(Touching_List), function(from) {
  14. if(is.null(Touching_List[[from]]==FALSE){      
  15. lines <- rgeos::gIntersection(Shapefile[from,], Shapefile[Touching_List[[from]],], byid = TRUE)
  16.       if (class(lines) == "SpatialCollections") {
  17.         list.Lines <- lapply(1:length(Touching_List[[from]]), function(to) {
  18.           line.single <- rgeos::gIntersection(Shapefile[from,], Shapefile[Touching_List[[from]][to],])
  19.           if (class(line.single) == "SpatialPoints") {
  20.             # Double the point to create a line
  21.             L1 <- rbind(line.single@coords, line.single@coords)
  22.             rownames(L1) <- letters[1:2]
  23.             Sl1 <- Line(L1)
  24.             Lines.single <- Lines(list(Sl1), ID = as.character(to))
  25.           } else if (class(line.single) == "SpatialLines") {
  26.             Lines.single <- line.single@lines[[1]]
  27.             Lines.single@ID <- as.character(to)
  28.           }
  29.           Lines.single
  30.         })
  31.         lines <- SpatialLines(list.Lines)
  32.       }
  33.       l_lines <- sp::SpatialLinesLengths(lines)
  34.       res <- data.frame(origin = from,
  35.                         perimeter = perimeters[from],
  36.                         touching = Touching_List[[from]],
  37.                         t.length = l_lines,
  38.                         t.pc = 100*l_lines/perimeters[from])
  39.       res
  40.     }
  41. })
  42.  
  43.     all.length.df <- do.call("rbind", all.length.list)
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy. OK, I Understand
 
Top