geometries

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#"
)
library(geometries)
library(Rcpp)

Header

#include "geometries/geometries.hpp"

Functions

SEXP make_geometries()

data as geometries

When one thinks of geometries in R, one of the most common data structures is the matrix. For example, in the sf world, an POINT is a single-row matrix (i.e, a vector)

sf::st_point( 1:2 )
POINT (1 2)

A LINESTRING is a matrix

sf::st_linestring( matrix( c(1,1,1,2,2,2,2,1,1,1), ncol = 2, byrow = T) )
LINESTRING (1 1, 1 2, 2 2, 2 1, 1 1)

and a POLYGON is a list of matrices

sf::st_polygon( list( matrix( c(1,1,1,2,2,2,2,1,1,1), ncol = 2, byrow = T) ) )
POLYGON ((1 1, 1 2, 2 2, 2 1, 1 1))

And to group these into a collection you would put each geometry inside a list

sf::st_sfc(
  list(
    sf::st_linestring( matrix( c(1,1,1,2,2,2,2,1,1,1), ncol = 2, byrow = T) )
    , sf::st_polygon( list( matrix( c(1,1,1,2,2,2,2,1,1,1), ncol = 2, byrow = T) ) )
  )
)

Geometry set for 2 features 
geometry type:  GEOMETRY
dimension:      XY
bbox:           xmin: 1 ymin: 1 xmax: 2 ymax: 2
CRS:            NA
LINESTRING (1 1, 1 2, 2 2, 2 1, 1 1)
POLYGON ((1 1, 1 2, 2 2, 2 1, 1 1))

From my limited research (i.e., practically none), I estimate most users will have a data.frame and will want to convert it into a collection of geometries.

For example, take a data.frame of x and y coordinates, and two id columns.

df <- data.frame(
  id1 = c( rep(1,12), rep(2, 12) )
  , id2 = c( rep(1:4, each = 3), rep(1:4, each = 3) )
  , x = 1:24
  , y = 24:1
)

df

You can think of the ID columns in this example as

Calling geometries::make_geometries() will split this data.frame by these ID columns and put the resulting matrices inside list elements.

cppFunction(

  depends = "geometries"
  , includes = '#include "geometries/geometries.hpp"'
  , code = '
    SEXP my_shape( SEXP df, SEXP id_cols, SEXP geometry_cols ) {
      return geometries::make_geometries( df, id_cols, geometry_cols );
    }
  '
  , plugins = "cpp11"
)

my_shape( df, c(0L,1L), c(2L,3L) )

# [[1]]
# [[1]][[1]]
#      [,1] [,2]
# [1,]    1   24
# [2,]    2   23
# [3,]    3   22
# 
# [[1]][[2]]
#      [,1] [,2]
# [1,]    4   21
# [2,]    5   20
# [3,]    6   19
# 
# [[1]][[3]]
#      [,1] [,2]
# [1,]    7   18
# [2,]    8   17
# [3,]    9   16
# 
# [[1]][[4]]
#      [,1] [,2]
# [1,]   10   15
# [2,]   11   14
# [3,]   12   13
# 
# 
# [[2]]
# [[2]][[1]]
#      [,1] [,2]
# [1,]   13   12
# [2,]   14   11
# [3,]   15   10
# 
# [[2]][[2]]
#      [,1] [,2]
# [1,]   16    9
# [2,]   17    8
# [3,]   18    7
# 
# [[2]][[3]]
#      [,1] [,2]
# [1,]   19    6
# [2,]   20    5
# [3,]   21    4
# 
# [[2]][[4]]
#      [,1] [,2]
# [1,]   22    3
# [2,]   23    2
# [3,]   24    1

Notice here there are no class attributes on the shapes. In {geometries} I only want to provide the tools to build these structures, then each user can define what they mean.

For example, if you want to define a class for each geometry you can supply a list containing a "class" vector as the 4th argument

cppFunction(
  depends = "geometries"
  , includes = '#include "geometries/geometries.hpp"'
  , code = '
    SEXP my_shape( Rcpp::DataFrame df, Rcpp::IntegerVector id_cols, Rcpp::IntegerVector geometry_cols, Rcpp::List class_attributes ) {
      return geometries::make_geometries( df, id_cols, geometry_cols, class_attributes );
    }
  '
  , plugins = "cpp11"
)

my_shape( df, c(0,1), c(2,3), list(class = "my_polygon") )

# [[1]]
# [[1]]
#      [,1] [,2]
# [1,]    1   24
# [2,]    2   23
# [3,]    3   22
# 
# [[2]]
#      [,1] [,2]
# [1,]    4   21
# [2,]    5   20
# [3,]    6   19
# 
# [[3]]
#      [,1] [,2]
# [1,]    7   18
# [2,]    8   17
# [3,]    9   16
# 
# [[4]]
#      [,1] [,2]
# [1,]   10   15
# [2,]   11   14
# [3,]   12   13
# 
# attr(,"class")
# [1] "my_polygon"
# 
# [[2]]
# [[1]]
#      [,1] [,2]
# [1,]   13   12
# [2,]   14   11
# [3,]   15   10
# 
# [[2]]
#      [,1] [,2]
# [1,]   16    9
# [2,]   17    8
# [3,]   18    7
# 
# [[3]]
#      [,1] [,2]
# [1,]   19    6
# [2,]   20    5
# [3,]   21    4
# 
# [[4]]
#      [,1] [,2]
# [1,]   22    3
# [2,]   23    2
# [3,]   24    1
# 
# attr(,"class")
# [1] "my_polygon"

Notice here that each list element now has a "my_polygon" class.

And if you have library(sf) loaded, setting the class as sfg POLYGON, you should see each element printed in the usual sf way

library(sf)
my_shape( df, c(0,1), c(2,3), list( class = c("XY", "POLYGON","sfg") ) )

# [[1]]
# POLYGON ((1 24, 2 23, 3 22), (4 21, 5 20, 6 19), (7 18, 8 17, 9 16), (10 15, 11 14, 12 13))
# 
# [[2]]
# POLYGON ((13 12, 14 11, 15 10), (16 9, 17 8, 18 7), (19 6, 20 5, 21 4), (22 3, 23 2, 24 1))

You can use this function to define any shape you want. The number of id columns you supply will determine how deeply nested the matrices are. If I add two more id columns, this will nest each matrix 2-levels deeper

df$id0 <- 1
df$id00 <- 1
head( df ) 

my_shape( df, c(0,1,4,5), c(2,3), list(class = "my_new_shape") )

# [[1]]
# [[1]]
# [[1]][[1]]
# [[1]][[1]][[1]]
#      [,1] [,2]
# [1,]    1   24
# [2,]    2   23
# [3,]    3   22
# 
# 
# 
# [[2]]
# [[2]][[1]]
# [[2]][[1]][[1]]
#      [,1] [,2]
# [1,]    4   21
# [2,]    5   20
# [3,]    6   19
# 
# 
# 
# [[3]]
# [[3]][[1]]
# [[3]][[1]][[1]]
#      [,1] [,2]
# [1,]    7   18
# [2,]    8   17
# [3,]    9   16
# 
# 
# 
# [[4]]
# [[4]][[1]]
# [[4]][[1]][[1]]
#      [,1] [,2]
# [1,]   10   15
# [2,]   11   14
# [3,]   12   13
# 
# 
# 
# attr(,"class")
# [1] "my_new_shape"
# 
# [[2]]
# [[1]]
# [[1]][[1]]
# [[1]][[1]][[1]]
#      [,1] [,2]
# [1,]   13   12
# [2,]   14   11
# [3,]   15   10
# 
# 
# 
# [[2]]
# [[2]][[1]]
# [[2]][[1]][[1]]
#      [,1] [,2]
# [1,]   16    9
# [2,]   17    8
# [3,]   18    7
# 
# 
# 
# [[3]]
# [[3]][[1]]
# [[3]][[1]][[1]]
#      [,1] [,2]
# [1,]   19    6
# [2,]   20    5
# [3,]   21    4
# 
# 
# 
# [[4]]
# [[4]][[1]]
# [[4]][[1]][[1]]
#      [,1] [,2]
# [1,]   22    3
# [2,]   23    2
# [3,]   24    1
# 
# 
# 
# attr(,"class")
# [1] "my_new_shape"


Try the geometries package in your browser

Any scripts or data that you put into this service are public.

geometries documentation built on Jan. 13, 2021, 8:12 p.m.