In ropenscilabs/geoops: 'GeoJSON' Topology Calculations and Operations

geoops

library("knitr")
hook_output <- knitr::knit_hooks$get("output")
knitr::knit_hooks$set(output = function(x, options) {
  lines <- options$output.lines
  if (is.null(lines)) {
    return(hook_output(x, options)) # pass to default hook
  }
  x <- unlist(strsplit(x, "\n"))
  more <- "..."
  if (length(lines) == 1) { # first n lines
    if (length(x) > lines) {
      # truncate the output, but add ....
      x <- c(head(x, lines), more)
    }
  } else {
    x <- c(
      if (abs(lines[1]) > 1) more else NULL,
      x[lines],
      if (length(x) > lines[abs(length(lines))]) more else NULL
    )
  }
  # paste these lines together
  x <- paste(c(x, ""), collapse = "\n")
  hook_output(x, options)
})

knitr::opts_chunk$set(
  comment = "#>",
  collapse = TRUE,
  warning = FALSE,
  message = FALSE,
  fig.path = "man/figures/"
)

geoops does spatial operations on GeoJSON.

geoops is inspired by the JS library turf (http://turfjs.org/). It's tagline is Advanced geospatial analysis for browsers and node. Turf works only with GeoJSON, as does geoops. I don't know JS that well, but it's easy enough to understand the language, so I've been porting Turf to C++ wrapped up in R. The C++ so we can have fast performance. We've also wrapped the Turf JS library itself in the package lawn (https://github.com/ropensci/lawn), but we should be able to get better performance out of C++.

geoops has a ways to go to include all the methods that Turf has, but we'll get there eventually.

All data is expected to be in WGS-84.

We use a library from Niels Lohmann (https://github.com/nlohmann/json) for working with JSON in C++.

See also:

• geojson: https://github.com/ropensci/geojson

Package API:

cat(paste(" -", paste(getNamespaceExports("geoops"), collapse = "\n - ")))

Installation

Stable version

install.packages("geoops")

Dev version

remotes::install_github("sckott/geoops")

library("geoops")

See the vignette (link here) to get started.

comparison to rgeos

FIXME!!! remove rgeos stuff as that pkg is gone

distance

pt1 <- '{"type":"Feature","properties":{"marker-color":"#f00"},"geometry":{"type":"Point","coordinates":[-75.343,39.984]}}'
pt2 <- '{"type":"Feature","properties":{"marker-color":"#0f0"},"geometry":{"type":"Point","coordinates":[-75.534,39.123]}}'
library(rgeos)
rgeospt1 <- rgeos::readWKT("POINT(0.5 0.5)")
rgeospt2 <- rgeos::readWKT("POINT(2 2)")

microbenchmark::microbenchmark(
  rgeos = rgeos::gDistance(rgeospt1, rgeospt2),
  geoops = geoops::geo_distance(pt1, pt2, units = "miles"),
  times = 10000L
)

nearest

point1 <- '{"type":["Feature"],"properties":{"marker-color":["#0f0"]},"geometry":{"type":["Point"],"coordinates":[28.9658,41.0101]}}'
point2 <- '{"type":["FeatureCollection"],"features":[{"type":"Feature","properties":{},"geometry":{"type":"Point","coordinates":[28.9739,41.0111]}},{"type":"Feature","properties":{},"geometry":{"type":"Point","coordinates":[28.9485,41.0242]}},{"type":"Feature","properties":{},"geometry":{"type":"Point","coordinates":[28.9387,41.0133]}}]}'
points <- '{"type":"FeatureCollection","features":[{"type":"Feature","properties":{},"geometry":{"type":"Point","coordinates":[28.9739,41.0111]}},{"type":"Feature","properties":{},"geometry":{"type":"Point","coordinates":[28.9485,41.0242]}},{"type":"Feature","properties":{},"geometry":{"type":"Point","coordinates":[28.9387,41.0133]}}]}'
g1 <- readWKT("MULTILINESTRING((34 54, 60 34), (0 10, 50 10, 100 50))")
g2 <- readWKT("POINT(100 30)")

microbenchmark::microbenchmark(
  rgeos = rgeos::gNearestPoints(g1, g2),
  geoops = geoops::geo_nearest(point1, points),
  times = 10000L
)

Example use case

expand

Get some GeoJSON data, a FeatureCollection of Polygons

file <- system.file("examples/zillow_or.geojson", package = "geoops")
x <- paste0(readLines(file), collapse = "")

Break each polygon into separate JSON string

library("jqr")
polys <- unclass(jq(x, ".features[]"))

Using geo_area, calculate the area of the polygon

areas <- vapply(polys, geo_area, 1, USE.NAMES = FALSE)

Visualize area of the polygons as a histogram

hist(areas, main = "")

Visualize some of the polygons, all of them

library(leaflet)
leaflet() %>%
  addProviderTiles(provider = "OpenStreetMap.Mapnik") %>%
  addGeoJSON(geojson = x) %>%
  setView(lng = -123, lat = 45, zoom = 7)

Just one of them

leaflet() %>%
  addProviderTiles(provider = "OpenStreetMap.Mapnik") %>%
  addGeoJSON(geojson = polys[1]) %>%
  setView(lng = -122.7, lat = 45.48, zoom = 13)

Meta

• Please report any issues or bugs.
• License: MIT
• Get citation information for geoops in R doing citation(package = 'geoops')
• Please note that this project is released with a Contributor Code of Conduct. By participating in this project you agree to abide by its terms.

ropenscilabs/geoops documentation built on Nov. 8, 2024, 4:01 a.m.