R version 3.6.2 (2019-12-12) -- "Dark and Stormy Night"
Copyright (C) 2019 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> suppressPackageStartupMessages(library(sf))
> # nc = st_read(system.file("gpkg/nc.gpkg", package="sf"))
> nc = st_read(system.file("shape/nc.shp", package="sf"), quiet = TRUE)
> nc_checked = st_transform(nc, 32119, check = TRUE)
> ncm = st_transform(nc, 32119)
> 
> x = st_transform(nc[1:10,], 32119)
> st_distance(x)
Units: [m]
           [,1]      [,2]      [,3]      [,4]      [,5]      [,6]     [,7]
 [1,]      0.00      0.00  25652.53 440561.15 299772.18 361529.53 419671.2
 [2,]      0.00      0.00      0.00 409428.54 268944.31 332589.77 388544.5
 [3,]  25652.53      0.00      0.00 367555.48 227017.51 290297.10 346667.9
 [4,] 440561.15 409428.54 367555.48      0.00  67226.86  45537.55      0.0
 [5,] 299772.18 268944.31 227017.51  67226.86      0.00      0.00  46527.4
 [6,] 361529.53 332589.77 290297.10  45537.55      0.00      0.00  30213.0
 [7,] 419671.16 388544.52 346667.93      0.00  46527.40  30213.00      0.0
 [8,] 384592.90 354294.11 312350.68  16130.34  11926.81      0.00      0.0
 [9,] 262354.29 231217.49 189310.35 140455.87      0.00  64606.24 119563.7
[10,]  71139.29  41943.83      0.00 330751.39 190182.39 252372.06 309862.0
           [,8]      [,9]     [,10]
 [1,] 384592.90 262354.29  71139.29
 [2,] 354294.11 231217.49  41943.83
 [3,] 312350.68 189310.35      0.00
 [4,]  16130.34 140455.87 330751.39
 [5,]  11926.81      0.00 190182.39
 [6,]      0.00  64606.24 252372.06
 [7,]      0.00 119563.73 309861.97
 [8,]      0.00  85533.12 275389.78
 [9,]  85533.12      0.00 152488.92
[10,] 275389.78 152488.92      0.00
> 
> st_is_valid(nc)
  [1] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
 [16] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
 [31] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
 [46] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
 [61] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
 [76] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
 [91] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
> 
> st_is_empty(st_sfc(st_point(), st_linestring()))
[1] TRUE TRUE
> 
> ops = c("intersects", #"disjoint", 
+ "touches", "crosses", "within", "contains", "overlaps", "equals", "covers", "covered_by", "equals_exact")
> for (op in ops) {
+ 	x = sf:::st_geos_binop(op, ncm[1:50,], ncm[51:100,], 0, NA_character_, FALSE)
+ 	x = sf:::st_geos_binop(op, ncm[1:50,], ncm[51:100,], 0, NA_character_, TRUE)
+ }
> 
> ops = c("intersects", #"disjoint", 
+ "touches", "crosses", "within", "contains", "overlaps", "covers", "covered_by")
> df = data.frame(ops = ops)
> df$equal = NA
> for (op in ops)
+ 	df[df$ops == op, "equal"] = identical(
+ 		sf:::st_geos_binop(op, ncm[1:50,], ncm[51:100,], 0, NA_character_, TRUE, FALSE),
+ 		sf:::st_geos_binop(op, ncm[1:50,], ncm[51:100,], 0, NA_character_, TRUE,  TRUE)
+ 	)
> df
         ops equal
1 intersects  TRUE
2    touches  TRUE
3    crosses  TRUE
4     within  TRUE
5   contains  TRUE
6   overlaps  TRUE
7     covers  TRUE
8 covered_by  TRUE
> 
> st_contains_properly(ncm[1:3,], ncm[1:3])
Sparse geometry binary predicate list of length 3, where the predicate was `contains_properly'
 1: (empty)
 2: (empty)
 3: (empty)
> 
> st_combine(nc)
Geometry set for 1 feature 
geometry type:  MULTIPOLYGON
dimension:      XY
bbox:           xmin: -84.32385 ymin: 33.88199 xmax: -75.45698 ymax: 36.58965
geographic CRS: NAD27
MULTIPOLYGON (((-81.47276 36.23436, -81.54084 3...
> 
> st_dimension(st_sfc(st_point(0:1), st_linestring(rbind(c(0,0),c(1,1))), 
+ 	st_polygon(list(rbind(c(0,0), c(1,0), c(1,1), c(0,1), c(0,0))))))
[1] 0 1 2
> 
> ncbb = st_as_sfc(st_bbox(nc))
> g = st_make_grid(ncbb)
although coordinates are longitude/latitude, st_relate_pattern assumes that they are planar
> x = st_intersection(nc, g)
although coordinates are longitude/latitude, st_intersection assumes that they are planar
Warning message:
attribute variables are assumed to be spatially constant throughout all geometries 
> x = st_intersection(g, nc)
although coordinates are longitude/latitude, st_intersection assumes that they are planar
> 
> ls = st_sfc(st_linestring(rbind(c(0,0),c(0,1))),
+ st_linestring(rbind(c(0,0),c(10,0))))
> 
> suppressWarnings(RNGversion("3.5.3"))
> set.seed(13531)
> 
> st_line_sample(ls, density = 1, type = "random")
Geometry set for 2 features 
geometry type:  MULTIPOINT
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 6.880179 ymax: 0.8878369
CRS:            NA
MULTIPOINT ((0 0.8878369))
MULTIPOINT ((0.2986488 0), (2.48417 0), (2.5678...
> 
> g = st_make_grid(ncbb, n = c(20,10))
although coordinates are longitude/latitude, st_relate_pattern assumes that they are planar
> 
> a1 = st_interpolate_aw(nc["BIR74"], g, FALSE)
although coordinates are longitude/latitude, st_intersection assumes that they are planar
Warning message:
In st_interpolate_aw.sf(nc["BIR74"], g, FALSE) :
  st_interpolate_aw assumes attributes are constant over areas of x
> sum(a1$BIR74) / sum(nc$BIR74) # not close to one: property is assumed spatially intensive
[1] 1.436169
> a2 = st_interpolate_aw(nc["BIR74"], g, extensive = TRUE)
although coordinates are longitude/latitude, st_intersection assumes that they are planar
Warning message:
In st_interpolate_aw.sf(nc["BIR74"], g, extensive = TRUE) :
  st_interpolate_aw assumes attributes are constant over areas of x
> sum(a2$BIR74) / sum(nc$BIR74)
[1] 1.000013
> 
> # missing x:
> g = st_make_grid(offset = c(0,0), cellsize = c(1,1), n = c(10,10))
> g = st_make_grid(what = "centers")
> length(g)
[1] 648
> g = st_make_grid(what = "corners")
> length(g)
[1] 703
> 
> g1 = st_make_grid(ncbb, 0.1, what = "polygons", square = FALSE)
although coordinates are longitude/latitude, st_relate_pattern assumes that they are planar
> g2 = st_make_grid(ncbb, 0.1, what = "points", square = FALSE)
although coordinates are longitude/latitude, st_intersects assumes that they are planar
although coordinates are longitude/latitude, st_intersects assumes that they are planar
> 
> # st_line_merge:
> mls = st_multilinestring(list(rbind(c(0,0), c(1,1)), rbind(c(2,0), c(1,1))))
> st_line_merge(mls)
LINESTRING (0 0, 1 1, 2 0)
> 
> if (sf_extSoftVersion()["GEOS"] >= "3.5.0") {
+  # voronoi:
+  set.seed(1)
+  x = st_multipoint(matrix(runif(10),,2))
+  box = st_polygon(list(rbind(c(0,0),c(1,0),c(1,1),c(0,1),c(0,0))))
+  v = st_sfc(st_voronoi(x, st_sfc(box)))
+  plot(v, col = 0, border = 1, axes = TRUE)
+  plot(box, add = TRUE, col = 0, border = 1) # a larger box is returned, as documented
+  plot(x, add = TRUE, col = 'red', cex=2, pch=16)
+  plot(st_intersection(st_cast(v), box)) # clip to smaller box
+  plot(x, add = TRUE, col = 'red', cex=2, pch=16)
+ 
+  v = st_voronoi(x)
+  print(class(v))
+  v = st_sfc(st_voronoi(st_sfc(x)))
+  print(class(v))
+  v = st_voronoi(st_sf(a = 1, geom = st_sfc(x)))
+  print(class(v))
+ }
[1] "XY"                 "GEOMETRYCOLLECTION" "sfg"               
[1] "sfc_GEOMETRYCOLLECTION" "sfc"                   
[1] "sf"         "data.frame"
> 
> i = st_intersects(ncm, ncm[1:88,])
> all.equal(i, t(t(i)))
[1] TRUE
> 
> # check use of pattern in st_relate:
> sfc = st_as_sfc(st_bbox(st_sfc(st_point(c(0,0)), st_point(c(3,3)))))
> grd = st_make_grid(sfc, n = c(3,3))
> st_intersects(grd)
Sparse geometry binary predicate list of length 9, where the predicate was `intersects'
 1: 1, 2, 4, 5
 2: 1, 2, 3, 4, 5, 6
 3: 2, 3, 5, 6
 4: 1, 2, 4, 5, 7, 8
 5: 1, 2, 3, 4, 5, 6, 7, 8, 9
 6: 2, 3, 5, 6, 8, 9
 7: 4, 5, 7, 8
 8: 4, 5, 6, 7, 8, 9
 9: 5, 6, 8, 9
> st_relate(grd, pattern = "****1****")
Sparse geometry binary predicate list of length 9, where the predicate was `relate_pattern'
 1: 1, 2, 4
 2: 1, 2, 3, 5
 3: 2, 3, 6
 4: 1, 4, 5, 7
 5: 2, 4, 5, 6, 8
 6: 3, 5, 6, 9
 7: 4, 7, 8
 8: 5, 7, 8, 9
 9: 6, 8, 9
> st_relate(grd, pattern = "****0****")
Sparse geometry binary predicate list of length 9, where the predicate was `relate_pattern'
 1: 5
 2: 4, 6
 3: 5
 4: 2, 8
 5: 1, 3, 7, 9
 6: 2, 8
 7: 5
 8: 4, 6
 9: 5
> st_rook = function(a, b = a, ...) st_relate(a, b, pattern = "F***1****", ...)
> st_rook(grd, sparse = FALSE)
       [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9]
 [1,] FALSE  TRUE FALSE  TRUE FALSE FALSE FALSE FALSE FALSE
 [2,]  TRUE FALSE  TRUE FALSE  TRUE FALSE FALSE FALSE FALSE
 [3,] FALSE  TRUE FALSE FALSE FALSE  TRUE FALSE FALSE FALSE
 [4,]  TRUE FALSE FALSE FALSE  TRUE FALSE  TRUE FALSE FALSE
 [5,] FALSE  TRUE FALSE  TRUE FALSE  TRUE FALSE  TRUE FALSE
 [6,] FALSE FALSE  TRUE FALSE  TRUE FALSE FALSE FALSE  TRUE
 [7,] FALSE FALSE FALSE  TRUE FALSE FALSE FALSE  TRUE FALSE
 [8,] FALSE FALSE FALSE FALSE  TRUE FALSE  TRUE FALSE  TRUE
 [9,] FALSE FALSE FALSE FALSE FALSE  TRUE FALSE  TRUE FALSE
> 
> #if (Sys.getenv("USER") %in% c("edzer", "travis")) { # memory leaks:
>   try(st_relate(st_point(), st_point(), pattern = "FF*FF****")) # error: use st_disjoint
Error in CPL_geos_binop(st_geometry(x), st_geometry(y), op, par, pattern,  : 
  use st_disjoint for this pattern
> #}
> 
> a = st_is_within_distance(nc[c(1:3,20),], nc[1:3,], 100000, sparse = FALSE)
> b = st_is_within_distance(nc[c(1:3,20),], nc[1:3,], units::set_units(100000, m), sparse = FALSE)
> all.equal(a, b)
[1] TRUE
> x = st_is_within_distance(nc[1:3,], nc[1:5,], 100000)
> y = st_is_within_distance(nc[1:3,], nc[1:5,], units::set_units(100, km))
> all.equal(x, y)
[1] TRUE
> 
> nc_3857 = st_transform(nc, 3857)
> a = st_is_within_distance(nc_3857[c(1:3,20),], nc_3857[1:3,], 100000, sparse = FALSE)
> b = st_is_within_distance(nc_3857[c(1:3,20),], nc_3857[1:3,], units::set_units(100000, m), sparse = FALSE)
> all.equal(a, b)
[1] TRUE
> x = st_is_within_distance(nc_3857, nc_3857, 100000)
> y = st_is_within_distance(nc_3857, nc_3857, units::set_units(100, km))
> all.equal(x, y)
[1] TRUE
> 
> pe = st_sfc(st_point())
> p = st_sfc(st_point(c(0,0)), st_point(c(0,1)), st_point(c(0,2)))
> st_distance(p, p)
     [,1] [,2] [,3]
[1,]    0    1    2
[2,]    1    0    1
[3,]    2    1    0
> st_distance(p, pe)
     [,1]
[1,]   NA
[2,]   NA
[3,]   NA
> st_distance(p, p, by_element = TRUE)
[1] 0 0 0
> st_crs(p) = 4326
> st_distance(p, p[c(2,3,1)], by_element = TRUE)
Units: [m]
[1] 110574.4 110575.1 221149.5
> p = st_transform(p, 3587)
> st_distance(p, p[c(2,3,1)], by_element = TRUE)
Units: [m]
[1] 144589.5 142873.3 287462.8
> 
> # from https://github.com/r-spatial/sf/issues/458 :
> pts <- st_sfc(st_point(c(.5,.5)), st_point(c(1.5, 1.5)), st_point(c(2.5, 2.5)))
> pol <- st_polygon(list(rbind(c(0,0), c(2,0), c(2,2), c(0,2), c(0,0))))
> pol_df <- data.frame(id = 1) 
> st_geometry(pol_df) <- st_sfc(pol)
> st_intersects(pts, pol_df[pol_df$id == 2,]) # with empty sf/sfc
Sparse geometry binary predicate list of length 3, where the predicate was `intersects'
 1: (empty)
 2: (empty)
 3: (empty)
> st_intersects(pts, pol_df[pol_df$id == 2,], sparse = FALSE) # with empty sf/sfc
    
[1,]
[2,]
[3,]
> 
> # st_node
> l = st_linestring(rbind(c(0,0), c(1,1), c(0,1), c(1,0), c(0,0)))
> st_node(l)
MULTILINESTRING ((0 0, 0.5 0.5), (0.5 0.5, 1 1, 0 1, 0.5 0.5), (0.5 0.5, 1 0, 0 0))
> st_node(st_sfc(l))
Geometry set for 1 feature 
geometry type:  MULTILINESTRING
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
MULTILINESTRING ((0 0, 0.5 0.5), (0.5 0.5, 1 1,...
> st_node(st_sf(a = 1, st_sfc(l)))
Simple feature collection with 1 feature and 1 field
geometry type:  MULTILINESTRING
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
  a                      st_sfc.l.
1 1 MULTILINESTRING ((0 0, 0.5 ...
> 
> # print.sgbp:
> (lst = st_disjoint(nc, nc))
although coordinates are longitude/latitude, st_intersects assumes that they are planar
Sparse geometry binary predicate list of length 100, where the predicate was `disjoint'
first 10 elements:
 1: 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, ...
 2: 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, ...
 3: 1, 4, 5, 6, 7, 8, 9, 11, 12, 13, ...
 4: 1, 2, 3, 5, 6, 8, 9, 10, 11, 12, ...
 5: 1, 2, 3, 4, 7, 8, 10, 11, 12, 13, ...
 6: 1, 2, 3, 4, 7, 9, 10, 11, 12, 13, ...
 7: 1, 2, 3, 5, 6, 9, 10, 11, 12, 13, ...
 8: 1, 2, 3, 4, 5, 9, 10, 11, 12, 13, ...
 9: 1, 2, 3, 4, 6, 7, 8, 10, 11, 12, ...
 10: 1, 2, 4, 5, 6, 7, 8, 9, 11, 13, ...
> # dim.sgbp:
> dim(lst)
[1] 100 100
> # as.matrix.sgbp:
> as.matrix(lst)[1:5, 1:5]
      [,1]  [,2]  [,3]  [,4]  [,5]
[1,] FALSE FALSE  TRUE  TRUE  TRUE
[2,] FALSE FALSE FALSE  TRUE  TRUE
[3,]  TRUE FALSE FALSE  TRUE  TRUE
[4,]  TRUE  TRUE  TRUE FALSE  TRUE
[5,]  TRUE  TRUE  TRUE  TRUE FALSE
> # negate:
> !lst
Sparse geometry binary predicate list of length 100, where the predicate was `!disjoint'
first 10 elements:
 1: 1, 2, 18, 19
 2: 1, 2, 3, 18
 3: 2, 3, 10, 18, 23, 25
 4: 4, 7, 56
 5: 5, 6, 9, 16, 28
 6: 5, 6, 8, 28
 7: 4, 7, 8, 17
 8: 6, 7, 8, 17, 20, 21
 9: 5, 9, 15, 16, 24, 31
 10: 3, 10, 12, 25, 26
> # as.data.frame:
> head(as.data.frame(lst), 10)
   row.id col.id
1       1      3
2       1      4
3       1      5
4       1      6
5       1      7
6       1      8
7       1      9
8       1     10
9       1     11
10      1     12
> 
> # snap:
> nc1 = st_transform(nc, 32119)
> g = st_make_grid(nc1, c(5000,5000), what = "centers")
> s = st_snap(nc1[1:3,], g, 2501*sqrt(2))
> sfg = st_snap(st_geometry(nc1)[[1]], g, 2501*sqrt(2))
> sfg = st_snap(st_geometry(nc1)[[1]], st_combine(g), 2501*sqrt(2))
> 
> # Hausdorff distance: http://geos.refractions.net/ro/doxygen_docs/html/classgeos_1_1algorithm_1_1distance_1_1DiscreteHausdorffDistance.html
> A = st_as_sfc("LINESTRING (0 0, 100 0, 10 100, 10 100)")
> B = st_as_sfc("LINESTRING (0 100, 0 10, 80 10)")
> st_distance(c(A,B))
         [,1]     [,2]
[1,] 0.000000 8.176236
[2,] 8.176236 0.000000
> st_distance(c(A,B), which = "Hausdorff")
         [,1]     [,2]
[1,]  0.00000 22.36068
[2,] 22.36068  0.00000
> st_distance(c(A,B), which = "Hausdorff", par = 0.001)
             [,1]         [,2]
[1,] 2.929643e-14 4.789000e+01
[2,] 4.789000e+01 2.131628e-14
> LE = st_as_sfc("LINESTRING EMPTY")
> st_distance(c(A, LE), which = "Hausdorff", par = 0.001)
             [,1] [,2]
[1,] 2.929643e-14   NA
[2,]           NA   NA
> 
> # one-argument st_intersection and st_difference:
> set.seed(131)
> m = rbind(c(0,0), c(1,0), c(1,1), c(0,1), c(0,0))
> p = st_polygon(list(m))
> n = 100
> l = vector("list", n)
> for (i in 1:n)
+    l[[i]] = p + 10 * runif(2)
> s = st_sfc(l)
> plot(s, col = sf.colors(categorical = TRUE, alpha = .5))
> d = st_difference(s) # sequential differences: s1, s2-s1, s3-s2-s1, ...
> plot(d, col = sf.colors(categorical = TRUE, alpha = .5))
> i = st_intersection(s) # all intersections
> plot(i, col = sf.colors(categorical = TRUE, alpha = .5))
> summary(lengths(st_overlaps(s, s)))
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
   0.00    2.00    3.50    3.66    5.00    8.00 
> summary(lengths(st_overlaps(d, d)))
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
      0       0       0       0       0       0 
> summary(lengths(st_overlaps(i, i)))
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
      0       0       0       0       0       0 
> 
> sf = st_sf(s)
> i = st_intersection(sf) # all intersections
> plot(i["n.overlaps"])
> summary(i$n.overlaps - lengths(i$origins))
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
      0       0       0       0       0       0 
> 
> # st_nearest_points:
> pt1 = st_point(c(.1,.1))
> pt2 = st_point(c(.9,.9))
> b1 = st_buffer(pt1, 0.1)
> b2 = st_buffer(pt2, 0.1)
> plot(b1, xlim = c(0,1), ylim = c(0,1))
> plot(b2, add = TRUE)
> (ls0 = try(st_nearest_points(b1, b2))) # sfg
Geometry set for 1 feature 
geometry type:  LINESTRING
dimension:      XY
bbox:           xmin: 0.1707107 ymin: 0.1707107 xmax: 0.8292893 ymax: 0.8292893
CRS:            NA
LINESTRING (0.1707107 0.1707107, 0.8292893 0.82...
> (ls = try(st_nearest_points(st_sfc(b1), st_sfc(b2)))) # sfc
Geometry set for 1 feature 
geometry type:  LINESTRING
dimension:      XY
bbox:           xmin: 0.1707107 ymin: 0.1707107 xmax: 0.8292893 ymax: 0.8292893
CRS:            NA
LINESTRING (0.1707107 0.1707107, 0.8292893 0.82...
> (ls = try(st_nearest_points(st_sfc(b1), st_sfc(b2), pairwise = TRUE))) # sfc
Geometry set for 1 feature 
geometry type:  LINESTRING
dimension:      XY
bbox:           xmin: 0.1707107 ymin: 0.1707107 xmax: 0.8292893 ymax: 0.8292893
CRS:            NA
LINESTRING (0.1707107 0.1707107, 0.8292893 0.82...
> identical(ls0, ls)
[1] TRUE
> # plot(ls, add = TRUE, col = 'red')
> 
> nc = read_sf(system.file("gpkg/nc.gpkg", package="sf"))
> plot(st_geometry(nc))
> ls = try(st_nearest_points(nc[1,], nc))
although coordinates are longitude/latitude, st_nearest_points assumes that they are planar
> # plot(ls, col = 'red', add = TRUE)
> pts = st_cast(ls, "POINT") # gives all start & end points
> # starting, "from" points, corresponding to x:
> plot(pts[seq(1, 200, 2)], add = TRUE, col = 'blue')
> # ending, "to" points, corresponding to y:
> plot(pts[seq(2, 200, 2)], add = TRUE, col = 'red')
> 
> # points to nearest features
> ls1 = st_linestring(rbind(c(0,0), c(1,0)))
> ls2 = st_linestring(rbind(c(0,0.1), c(1,0.1)))
> ls3 = st_linestring(rbind(c(0,1), c(1,1)))
> (l = st_sfc(ls1, ls2, ls3))
Geometry set for 3 features 
geometry type:  LINESTRING
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
LINESTRING (0 0, 1 0)
LINESTRING (0 0.1, 1 0.1)
LINESTRING (0 1, 1 1)
> 
> p1 = st_point(c(0.1, -0.1))
> p2 = st_point(c(0.1, 0.11))
> p3 = st_point(c(0.1, 0.09))
> p4 = st_point(c(0.1, 0.9))
> p5 = st_point()
> 
> (p = st_sfc(p1, p2, p3, p4, p5))
Geometry set for 5 features  (with 1 geometry empty)
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0.1 ymin: -0.1 xmax: 0.1 ymax: 0.9
CRS:            NA
POINT (0.1 -0.1)
POINT (0.1 0.11)
POINT (0.1 0.09)
POINT (0.1 0.9)
POINT EMPTY
> #st_nearest_points(p, l)
> n = try(st_nearest_feature(p,l))
> if (!inherits(n, "try-error")) {
+   print(st_nearest_points(p, l[n], pairwise = TRUE))
+   print(st_nearest_feature(p, l))
+   print(st_nearest_feature(p, st_sfc()))
+   print(st_nearest_feature(st_sfc(), l))
+ }
Geometry set for 5 features  (with 1 geometry empty)
geometry type:  LINESTRING
dimension:      XY
bbox:           xmin: 0.1 ymin: -0.1 xmax: 0.1 ymax: 1
CRS:            NA
LINESTRING (0.1 -0.1, 0.1 0)
LINESTRING (0.1 0.11, 0.1 0.1)
LINESTRING (0.1 0.09, 0.1 0.1)
LINESTRING (0.1 0.9, 0.1 1)
LINESTRING EMPTY
[1]  1  2  2  3 NA
[1] NA NA NA NA NA
integer(0)
> 
> # can do centroid of empty geom:
> st_centroid(st_polygon())
POINT EMPTY
> 
> #999:
> pt = data.frame(x=1:2, y=1:2,a=letters[1:2])
> pt = st_as_sf(pt, coords=c("x","y"))
> 
> bf =st_buffer(pt, dist=0.3)
> 
> st_within(pt,bf, sparse=FALSE)
      [,1]  [,2]
[1,]  TRUE FALSE
[2,] FALSE  TRUE
> st_within(pt[1,], bf[1,], sparse = FALSE)
     [,1]
[1,] TRUE
> st_relate(pt[1,], bf[1,], pattern = "T*F**F***", sparse = FALSE)
     [,1]
[1,] TRUE
> 
> sf:::is_symmetric(pattern = "010121010")
[1] TRUE
> sf:::is_symmetric(pattern = "010121021")
[1] FALSE
> 
> st_intersects(st_point(0:1), st_point(2:3)) # sfg method
Sparse geometry binary predicate list of length 1, where the predicate was `intersects'
 1: (empty)
> 
> if (sf_extSoftVersion()["GEOS"] >= "3.7.0") {
+ 	ls = st_linestring(rbind(c(1,1), c(2,2), c(3,3)))
+ 	print(st_reverse(ls))
+ 	print(st_reverse(st_sfc(ls)))
+ 	print(st_reverse(st_sf(a = 2, geom = st_sfc(ls))))
+ }
LINESTRING (3 3, 2 2, 1 1)
Geometry set for 1 feature 
geometry type:  LINESTRING
dimension:      XY
bbox:           xmin: 1 ymin: 1 xmax: 3 ymax: 3
CRS:            NA
LINESTRING (3 3, 2 2, 1 1)
Simple feature collection with 1 feature and 1 field
geometry type:  LINESTRING
dimension:      XY
bbox:           xmin: 1 ymin: 1 xmax: 3 ymax: 3
CRS:            NA
  a                       geom
1 2 LINESTRING (3 3, 2 2, 1 1)
> 
> proc.time()
   user  system elapsed 
  4.369   0.075   4.463