hypertidy / toposhop

2 stars 0 forks source link

maps - restore arc-node #6

Open mdsumner opened 7 years ago

mdsumner commented 7 years ago

The maps package creates its objects from arc-node topology, or something very like it i.e. TopoJSON without the differencing or offset/scaling. That's why the default output is not understandable in terms of groupings.

mdsumner commented 7 years ago

topojson-parse branch contains test code

mdsumner commented 7 years ago

## dig under map()
mp <- maps:::map.poly(database)
plot(mp[1:2], pch = ".")
xlim <- c(-1e+30, 1e+30)
ylim <- c(-1e+30, 1e+30)
## dig under map.poly()
gon <- maps:::mapname(database, patterns  = "", exact = FALSE)
mg  <- maps:::mapgetg(database, gon, fill = FALSE , xlim, ylim)
coord <- maps:::mapgetl(database, mg$number, xlim, ylim, fill = TRUE) 
plot(coord[1:2], pch = ".")

database <- "world" 
dbname <- paste(database, "MapEnv", sep = "")
lines <- maps:::mapname(database, patterns = ".", exact = FALSE)
mapbase <- paste(Sys.getenv(get(dbname)), database, sep = "")
as.polygon <- FALSE

line <- maps:::mapgetg(database, gons = lines, fill = as.polygon, xlim = xlim, ylim = ylim)
str(line)
# List of 3
# $ number: int [1:2366] 1 2 3 4 5 6 7 8 9 10 ...
# $ size  : int [1:1639] 1 6 4 3 1 5 1 1 1 2 ...
# $ name  : chr [1:1639] "Aruba" "Afghanistan" "Angola" "Angola:Cabinda" ...
sum(line$size)
# [2] 2366

nline <- as.integer(length(lines))
fill <- TRUE
z <- .C("mapgetl", PACKAGE = "maps", as.character(mapbase), 
        linesize = as.integer(lines), error = as.integer(nline), 
        as.integer(0), as.double(0), as.double(0), as.double(c(xlim, ylim)), 
        as.integer(fill))[c("linesize", "error")]                                                      
str(z)
# List of 2
# $ linesize: int [1:1639] 10 12 157 67 64 12 103 51 66 62 ...
# $ error   : int 1639

sum(z$linesize)
#[1] 66411

ok <- z$linesize != 0
lines <- lines[ok]
nline <- length(lines)
if (nline == 0) 
  return(integer(0))
linesize <- z$linesize[ok]
N <- sum(linesize) + nline - 1
coord <- .C("mapgetl", PACKAGE = "maps", as.character(mapbase), 
            as.integer(lines), as.integer(nline), as.integer(1), 
            x = double(N), y = double(N), range = double(4), 
            as.integer(fill))  ##[c("x", "y", "range")]
str(coord)
# List of 8
# $      : chr "C:/Users/mdsumner/Documents/R/win-library/3.3/maps/mapdata/world2"
# $      : int [1:1639] 1 2 3 4 5 6 7 8 9 10 ...
# $      : int 1639
# $      : int 1
# $ x    : num [1:68049] 290 290 290 290 290 ...
# $ y    : num [1:68049] 12.5 12.4 12.4 12.5 12.5 ...
# $ range: num [1:4] 0 360 -85.2 83.6
# $      : int 0
library(dplyr)
v <- tibble::tibble(x_ = coord$x, y_ = coord$y) %>% 
  mutate(vertex_ = row_number(), group = cumsum(is.na(x_)) + 1) %>% 
  filter(!is.na(x_)) %>% mutate(group = as.integer(factor(group)))

ggplot(v, aes(x = x_, y = y_, group = group)) + geom_point(pch = ".")

arc <- tibble(vertex_ = v$vertex_, arc = rep(seq_along(line$number), line$size))

makepoly takes coord, z$linesize
# gon <- mapname(database, regions, exact)
# n <- length(gon)
# if (n == 0) 
#   stop("no recognized region names")
# if (is.null(xlim)) 
#   xlim <- c(-1e+30, 1e+30)
# if (is.null(ylim)) 
#   ylim <- c(-1e+30, 1e+30)
# line <- mapgetg(database, gon, as.polygon, xlim, ylim)
# if (length(line$number) == 0) 
#   stop("nothing to draw: all regions out of bounds")
# if (as.polygon) {
#   coord <- mapgetl(database, line$number, xlim, ylim, 
#                    fill)
#   gonsize <- line$size
#   keep <- rep(TRUE, length(gonsize))
#   coord[c("x", "y")] <- makepoly(coord, gonsize, keep)
# }
# else {
#   l <- abs(line$number)
#   if (boundary && interior) 
#     l <- unique(l)
#   else if (boundary) 
#     l <- l[!match(l, l[duplicated(l)], FALSE)]
#   else l <- l[duplicated(l)]
#   coord <- mapgetl(database, l, xlim, ylim, fill)
#   if (length(coord) == 0) 
#     stop("all data out of bounds")
# }