R version 4.0.0 (2020-04-24) -- "Arbor Day"
Copyright (C) 2020 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))
> library(testthat)
> 
> p = st_point(c(1/3,1/6))
> st_sfc(p, precision = 1000)
Geometry set for 1 feature 
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0.3333333 ymin: 0.1666667 xmax: 0.3333333 ymax: 0.1666667
CRS:            NA
precision:      1000 
POINT (0.3333333 0.1666667)
> st_as_sfc(st_as_binary(st_sfc(p, precision = 0L)))
Geometry set for 1 feature 
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0.3333333 ymin: 0.1666667 xmax: 0.3333333 ymax: 0.1666667
CRS:            NA
POINT (0.3333333 0.1666667)
> st_as_sfc(st_as_binary(st_sfc(p, precision = 1000)))
Geometry set for 1 feature 
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0.333 ymin: 0.167 xmax: 0.333 ymax: 0.167
CRS:            NA
POINT (0.333 0.167)
> st_as_sfc(st_as_binary(st_sfc(p, precision = 1000000)))
Geometry set for 1 feature 
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0.333333 ymin: 0.166667 xmax: 0.333333 ymax: 0.166667
CRS:            NA
POINT (0.333333 0.166667)
> st_as_sfc(st_as_binary(st_sfc(p, precision = 10L)))
Geometry set for 1 feature 
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0.3 ymin: 0.2 xmax: 0.3 ymax: 0.2
CRS:            NA
POINT (0.3 0.2)
> st_as_sfc(st_as_binary(st_sfc(p, precision = -1)))
Geometry set for 1 feature 
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0.3333333 ymin: 0.1666667 xmax: 0.3333333 ymax: 0.1666667
CRS:            NA
POINT (0.3333333 0.1666667)
> 
> d = data.frame(a = 1:2)
> d$geom = c("POINT(0 0)", "POINT(1 1)")
> 
> st_as_sf(d, wkt = "geom")
Simple feature collection with 2 features and 1 field
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
  a        geom
1 1 POINT (0 0)
2 2 POINT (1 1)
> st_as_sf(d, wkt = 2)
Simple feature collection with 2 features and 1 field
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
  a        geom
1 1 POINT (0 0)
2 2 POINT (1 1)
> st_as_sf(d, wkt = "geom", remove = FALSE)
Simple feature collection with 2 features and 2 fields
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
  a       geom    geometry
1 1 POINT(0 0) POINT (0 0)
2 2 POINT(1 1) POINT (1 1)
> 
> st_as_sfc(c("POINT(0 0)", "POINT(1 1)"))
Geometry set for 2 features 
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
POINT (0 0)
POINT (1 1)
> st_as_sfc(c("POINT(0 0)", "POINT(1 1)", "POLYGON((0 0,1 1,0 1,0 0))"))
Geometry set for 3 features 
geometry type:  GEOMETRY
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
POINT (0 0)
POINT (1 1)
POLYGON ((0 0, 1 1, 0 1, 0 0))
> st_as_sfc(character(0))
Geometry set for 0 features 
bbox:           xmin: NA ymin: NA xmax: NA ymax: NA
CRS:            NA
> x = st_as_sfc(character(0), 4326)
> y = st_as_sfc(character(0), crs = 4326)
> all.equal(x, y)
[1] TRUE
> st_as_sfc(c("POINT(0 0)", "POINT(1 1)", "POLYGON((0 0,1 1,0 1,0 0))"),
+ 	"+proj=longlat +datum=WGS84")
Geometry set for 3 features 
geometry type:  GEOMETRY
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            +proj=longlat +datum=WGS84
POINT (0 0)
POINT (1 1)
POLYGON ((0 0, 1 1, 0 1, 0 0))
> dg = st_as_sf(d, wkt = "geom")
> print(dg, n = 1)
Simple feature collection with 2 features and 1 field
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
First 1 features:
  a        geom
1 1 POINT (0 0)
> head(st_as_sf(d, wkt = "geom"), 1)
Simple feature collection with 1 feature and 1 field
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 0 ymax: 0
CRS:            NA
  a        geom
1 1 POINT (0 0)
> 
> d$geom = st_as_sfc(d$geom)
> d1 = d
> attr(d1, "sf_col") = "geom"
> st_geometry(d1) = d$geom
> 
> d$geometry = d$geom # second geometry list-column
> expect_warning(st_geometry(d) <- d$geom)
> d
Simple feature collection with 2 features and 1 field
Active geometry column: geom
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
  a        geom    geometry
1 1 POINT (0 0) POINT (0 0)
2 2 POINT (1 1) POINT (1 1)
> 
> x = st_sfc(list(st_point(0:1), st_point(0:1)), crs = 4326)
> # don't warn when replacing crs with identical value:
> st_sfc(x, crs = 4326)
Geometry set for 2 features 
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 1 xmax: 0 ymax: 1
geographic CRS: WGS 84
POINT (0 1)
POINT (0 1)
> y = st_sfc(x, crs = "+proj=longlat +datum=WGS84 +no_defs")
> # but do when it changes:
> y = st_sfc(x, crs = 3857)
Warning message:
st_crs<- : replacing crs does not reproject data; use st_transform for that 
> 
> p = st_point(0:1)
> st_cast(p, "MULTIPOINT")
MULTIPOINT ((0 1))
> mp = st_multipoint(rbind(c(0,1), c(2,2)))
> st_cast(mp, "POINT")
POINT (0 1)
Warning message:
In st_cast.MULTIPOINT(mp, "POINT") : point from first coordinate only
> st_cast(mp, "MULTIPOINT")
MULTIPOINT ((0 1), (2 2))
> 
> # geometry collection to its elements:
> st_cast(st_geometrycollection(list(mp)), "POINT")
POINT (0 1)
Warning message:
In st_cast.MULTIPOINT(x[[1]], to, ...) : point from first coordinate only
> st_cast(st_geometrycollection(list(mp)), "MULTIPOINT")
MULTIPOINT ((0 1), (2 2))
> st_cast(st_geometrycollection(list(p,mp)), "MULTIPOINT")
MULTIPOINT ((0 1))
Warning message:
In st_cast.GEOMETRYCOLLECTION(st_geometrycollection(list(p, mp)),  :
  only first part of geometrycollection is retained
> 
> mp = st_multipoint(rbind(c(0,1)))
> x = st_sfc(p, mp)
> st_cast(x, "POINT")
Geometry set for 2 features 
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 1 xmax: 0 ymax: 1
CRS:            NA
POINT (0 1)
POINT (0 1)
Warning message:
In st_cast.MULTIPOINT(X[[i]], ...) : point from first coordinate only
> 
> sf = st_sf(a = 3:2, geom = x)
> st_cast(sf, "POINT")
Simple feature collection with 2 features and 1 field
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 1 xmax: 0 ymax: 1
CRS:            NA
  a        geom
1 3 POINT (0 1)
2 2 POINT (0 1)
Warning message:
In st_cast.MULTIPOINT(X[[i]], ...) : point from first coordinate only
> 
> suppressPackageStartupMessages( library(dplyr) )
> 
> x %>% st_cast("POINT")
Geometry set for 2 features 
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 1 xmax: 0 ymax: 1
CRS:            NA
POINT (0 1)
POINT (0 1)
Warning message:
In st_cast.MULTIPOINT(X[[i]], ...) : point from first coordinate only
> 
> # points:
> mp = st_multipoint(rbind(c(0,1))) # single-point multipoint
> st_sfc(p,mp) %>% st_cast("POINT")
Geometry set for 2 features 
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 1 xmax: 0 ymax: 1
CRS:            NA
POINT (0 1)
POINT (0 1)
Warning message:
In st_cast.MULTIPOINT(X[[i]], ...) : point from first coordinate only
> st_sfc(p,mp) %>% st_cast("MULTIPOINT")
Geometry set for 2 features 
geometry type:  MULTIPOINT
dimension:      XY
bbox:           xmin: 0 ymin: 1 xmax: 0 ymax: 1
CRS:            NA
MULTIPOINT ((0 1))
MULTIPOINT ((0 1))
> 
> # lines:
> pts = rbind(c(0,0), c(1,1), c(2,1))
> st_sfc(st_linestring(pts), st_multilinestring(list(pts))) %>% st_cast("LINESTRING")
Geometry set for 2 features 
geometry type:  LINESTRING
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 2 ymax: 1
CRS:            NA
LINESTRING (0 0, 1 1, 2 1)
LINESTRING (0 0, 1 1, 2 1)
> st_sfc(st_linestring(pts), st_multilinestring(list(pts))) %>% st_cast("MULTILINESTRING")
Geometry set for 2 features 
geometry type:  MULTILINESTRING
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 2 ymax: 1
CRS:            NA
MULTILINESTRING ((0 0, 1 1, 2 1))
MULTILINESTRING ((0 0, 1 1, 2 1))
> 
> # polygons:
> pts = rbind(c(0,0), c(1,1), c(0,1), c(0,0))
> st_sfc(st_polygon(list(pts)), st_multipolygon(list(list(pts)))) %>% st_cast("POLYGON")
Geometry set for 2 features 
geometry type:  POLYGON
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
POLYGON ((0 0, 1 1, 0 1, 0 0))
POLYGON ((0 0, 1 1, 0 1, 0 0))
> st_sfc(st_polygon(list(pts)), st_multipolygon(list(list(pts)))) %>% st_cast("MULTIPOLYGON")
Geometry set for 2 features 
geometry type:  MULTIPOLYGON
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
MULTIPOLYGON (((0 0, 1 1, 0 1, 0 0)))
MULTIPOLYGON (((0 0, 1 1, 0 1, 0 0)))
> 
> 
> st_sfc(st_geometrycollection(list(p)), st_geometrycollection(list(mp))) %>% st_cast()
Geometry set for 2 features 
geometry type:  GEOMETRY
dimension:      XY
bbox:           xmin: 0 ymin: 1 xmax: 0 ymax: 1
CRS:            NA
POINT (0 1)
MULTIPOINT ((0 1))
> st_sfc(st_geometrycollection(list(p)), st_geometrycollection(list(mp))) %>%
+ 	st_cast() %>%
+ 	st_cast("POINT")
Geometry set for 2 features 
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 1 xmax: 0 ymax: 1
CRS:            NA
POINT (0 1)
POINT (0 1)
Warning message:
In st_cast.MULTIPOINT(X[[i]], ...) : point from first coordinate only
> 
> p = rbind(c(0,0),c(1,0),c(1,1),c(0,1),c(0,0))
> pol = st_polygon(list(p))
> # plot(pol)
> try(plot(st_polygonize(pol))) # --> breaks
Error in st_polygonize.sfc(st_sfc(x)) : 
  inherits(x, "sfc_LINESTRING") || inherits(x, "sfc_MULTILINESTRING") is not TRUE
> st_length(st_sfc(st_point(c(0,0))))
[1] 0
> 
> try(as(st_sfc(st_linestring(matrix(1:9,3))), "Spatial"))
Error in StopZ(zm) : 
  sp supports Z dimension only for POINT and MULTIPOINT.
use `st_zm(...)` to coerce to XY dimensions
> 
> # check conus is present:
> x = st_sfc(st_point(c(-90,35)), st_point(c(-80,36)),
+ 	crs = "+proj=longlat +datum=NAD27")
> y = st_transform(x, 3857)
> 
> sf_extSoftVersion()[1:3]
   GEOS    GDAL  proj.4 
"3.8.0" "3.0.4" "7.0.0" 
> 
> # Ops.sfc:
> ls = st_sfc(st_linestring(rbind(c(0,0),c(0,1))))
> ls * 2
Geometry set for 1 feature 
geometry type:  LINESTRING
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 0 ymax: 2
CRS:            NA
LINESTRING (0 0, 0 2)
> ls - 2
Geometry set for 1 feature 
geometry type:  LINESTRING
dimension:      XY
bbox:           xmin: -2 ymin: -2 xmax: -2 ymax: -1
CRS:            NA
LINESTRING (-2 -2, -2 -1)
> (ls + 2) %% 3
Geometry set for 1 feature 
geometry type:  LINESTRING
dimension:      XY
bbox:           xmin: 2 ymin: 0 xmax: 2 ymax: 2
CRS:            NA
LINESTRING (2 2, 2 0)
> ls / ls
Geometry set for 1 feature  (with 1 geometry empty)
geometry type:  GEOMETRYCOLLECTION
dimension:      XY
bbox:           xmin: NA ymin: NA xmax: NA ymax: NA
CRS:            NA
GEOMETRYCOLLECTION EMPTY
> p_ = st_point(0:1)
> ll = st_sfc(ls[[1]], p_)
> ll & st_sfc(p_)
Geometry set for 2 features 
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 1 xmax: 0 ymax: 1
CRS:            NA
POINT (0 1)
POINT (0 1)
> ll | st_sfc(p_)
Geometry set for 2 features 
geometry type:  GEOMETRY
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 0 ymax: 1
CRS:            NA
LINESTRING (0 0, 0 1)
POINT (0 1)
> ll %/% st_sfc(p_)
Geometry set for 2 features  (with 1 geometry empty)
geometry type:  GEOMETRY
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 0 ymax: 1
CRS:            NA
LINESTRING (0 0, 0 1)
GEOMETRYCOLLECTION EMPTY
> ll == st_sfc(p_)
[1] FALSE  TRUE
> ll != st_sfc(p_)
[1]  TRUE FALSE
> 
> 
> str(x)
sfc_POINT of length 2; first list element:  'XY' num [1:2] -90 35
> nc = st_read(system.file("shape/nc.shp", package="sf"), quiet = TRUE)
> str(nc)
Classes 'sf' and 'data.frame':	100 obs. of  15 variables:
 $ AREA     : num  0.114 0.061 0.143 0.07 0.153 0.097 0.062 0.091 0.118 0.124 ...
 $ PERIMETER: num  1.44 1.23 1.63 2.97 2.21 ...
 $ CNTY_    : num  1825 1827 1828 1831 1832 ...
 $ CNTY_ID  : num  1825 1827 1828 1831 1832 ...
 $ NAME     : chr  "Ashe" "Alleghany" "Surry" "Currituck" ...
 $ FIPS     : chr  "37009" "37005" "37171" "37053" ...
 $ FIPSNO   : num  37009 37005 37171 37053 37131 ...
 $ CRESS_ID : int  5 3 86 27 66 46 15 37 93 85 ...
 $ BIR74    : num  1091 487 3188 508 1421 ...
 $ SID74    : num  1 0 5 1 9 7 0 0 4 1 ...
 $ NWBIR74  : num  10 10 208 123 1066 ...
 $ BIR79    : num  1364 542 3616 830 1606 ...
 $ SID79    : num  0 3 6 2 3 5 2 2 2 5 ...
 $ NWBIR79  : num  19 12 260 145 1197 ...
 $ geometry :sfc_MULTIPOLYGON of length 100; first list element: List of 1
  ..$ :List of 1
  .. ..$ : num [1:27, 1:2] -81.5 -81.5 -81.6 -81.6 -81.7 ...
  ..- attr(*, "class")= chr [1:3] "XY" "MULTIPOLYGON" "sfg"
 - attr(*, "sf_column")= chr "geometry"
 - attr(*, "agr")= Factor w/ 3 levels "constant","aggregate",..: NA NA NA NA NA NA NA NA NA NA ...
  ..- attr(*, "names")= chr [1:14] "AREA" "PERIMETER" "CNTY_" "CNTY_ID" ...
> bb = st_as_sfc(st_bbox(nc))
> format(st_bbox(nc))
[1] "((-84.32385,33.88199),(-75.45698,36.58965))"
> 
> st_agr("constant")
[1] constant
Levels: constant aggregate identity
> st_agr()
[1] <NA>
Levels: constant aggregate identity
> x <- st_sf(a = 1:2, b = 3:4, geom = x, agr = c("constant", "aggregate"))
> suppressPackageStartupMessages(library(dplyr))
> y <- x %>% st_set_agr("constant")
> y
Simple feature collection with 2 features and 2 fields
Attribute-geometry relationship: 2 constant, 0 aggregate, 0 identity
geometry type:  POINT
dimension:      XY
bbox:           xmin: -90 ymin: 35 xmax: -80 ymax: 36
CRS:            +proj=longlat +datum=NAD27
  a b           geom
1 1 3 POINT (-90 35)
2 2 4 POINT (-80 36)
> 
> sf1 <- st_sf(a = c("x", "y"), geom = st_sfc(st_point(3:4), st_point(3:4)))
> sf1[names(sf1)]
Simple feature collection with 2 features and 1 field
geometry type:  POINT
dimension:      XY
bbox:           xmin: 3 ymin: 4 xmax: 3 ymax: 4
CRS:            NA
  a        geom
1 x POINT (3 4)
2 y POINT (3 4)
> 
> st_bbox(sf1)
xmin ymin xmax ymax 
   3    4    3    4 
> bb = st_bbox(nc)
> bb
     xmin      ymin      xmax      ymax 
-84.32385  33.88199 -75.45698  36.58965 
> st_crs(bb)
Coordinate Reference System:
  User input: NAD27 
  wkt:
GEOGCRS["NAD27",
    DATUM["North American Datum 1927",
        ELLIPSOID["Clarke 1866",6378206.4,294.978698213898,
            LENGTHUNIT["metre",1]]],
    PRIMEM["Greenwich",0,
        ANGLEUNIT["degree",0.0174532925199433]],
    CS[ellipsoidal,2],
        AXIS["latitude",north,
            ORDER[1],
            ANGLEUNIT["degree",0.0174532925199433]],
        AXIS["longitude",east,
            ORDER[2],
            ANGLEUNIT["degree",0.0174532925199433]],
    ID["EPSG",4267]]
> st_bbox(c(xmin = 16.1, xmax = 16.6, ymin = 48.6, ymax = 47.9), crs = st_crs(4326))
xmin ymin xmax ymax 
16.1 48.6 16.6 47.9 
> st_bbox(c(xmin = 16.1, xmax = 16.6, ymin = 48.6, ymax = 47.9), crs = 4326)
xmin ymin xmax ymax 
16.1 48.6 16.6 47.9 
> 
> bb$xrange
     xmin      xmax 
-84.32385 -75.45698 
> bb$yrange
    ymin     ymax 
33.88199 36.58965 
> bb$xmin
     xmin 
-84.32385 
> bb$ymin
    ymin 
33.88199 
> bb$xmax
     xmax 
-75.45698 
> bb$ymax
    ymax 
36.58965 
> try(bb$foo)
Error in `$.bbox`(bb, foo) : unsupported name
> 
> # merge:
> a = data.frame(a = 1:3, b = 5:7)
> st_geometry(a) = st_sfc(st_point(c(0,0)), st_point(c(1,1)), st_point(c(2,2)))
> b = data.frame(x = c("a", "b", "c"), b = c(2,5,6))
> merge(a, b)
Simple feature collection with 2 features and 3 fields
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
  b a x    geometry
1 5 1 b POINT (0 0)
2 6 2 c POINT (1 1)
> merge(a, b, all = TRUE)
Simple feature collection with 4 features and 3 fields (with 1 geometry empty)
geometry type:  GEOMETRY
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 2 ymax: 2
CRS:            NA
  b  a    x                 geometry
1 2 NA    a GEOMETRYCOLLECTION EMPTY
2 5  1    b              POINT (0 0)
3 6  2    c              POINT (1 1)
4 7  3 <NA>              POINT (2 2)
> 
> # joins:
> inner_join(a, b)
Joining, by = "b"
Simple feature collection with 2 features and 3 fields
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
  a b x    geometry
1 1 5 b POINT (0 0)
2 2 6 c POINT (1 1)
> left_join(a, b)
Joining, by = "b"
Simple feature collection with 3 features and 3 fields
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 2 ymax: 2
CRS:            NA
  a b    x    geometry
1 1 5    b POINT (0 0)
2 2 6    c POINT (1 1)
3 3 7 <NA> POINT (2 2)
> right_join(a, b)
Joining, by = "b"
Simple feature collection with 3 features and 3 fields (with 1 geometry empty)
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
   a b x    geometry
1 NA 2 a POINT EMPTY
2  1 5 b POINT (0 0)
3  2 6 c POINT (1 1)
> full_join(a, b)
Joining, by = "b"
Simple feature collection with 4 features and 3 fields (with 1 geometry empty)
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 2 ymax: 2
CRS:            NA
   a b    x    geometry
1  1 5    b POINT (0 0)
2  2 6    c POINT (1 1)
3  3 7 <NA> POINT (2 2)
4 NA 2    a POINT EMPTY
> semi_join(a, b)
Joining, by = "b"
Simple feature collection with 2 features and 2 fields
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
  a b    geometry
1 1 5 POINT (0 0)
2 2 6 POINT (1 1)
> anti_join(a, b)
Joining, by = "b"
Simple feature collection with 1 feature and 2 fields
geometry type:  POINT
dimension:      XY
bbox:           xmin: 2 ymin: 2 xmax: 2 ymax: 2
CRS:            NA
  a b    geometry
1 3 7 POINT (2 2)
> left_join(a, data.frame(b, geometry = 1), by = "b")
Simple feature collection with 3 features and 4 fields
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 2 ymax: 2
CRS:            NA
  a b    x geometry.y  geometry.x
1 1 5    b          1 POINT (0 0)
2 2 6    c          1 POINT (1 1)
3 3 7 <NA>         NA POINT (2 2)
> 
> # st_joins:
> a = st_sf(a = 1:3,
+  geom = st_sfc(st_point(c(1,1)), st_point(c(2,2)), st_point(c(3,3))))
> b = st_sf(a = 11:14,
+  geom = st_sfc(st_point(c(10,10)), st_point(c(2,2)), st_point(c(2,2)), st_point(c(3,3))))
> st_join(a, b)
Simple feature collection with 4 features and 2 fields
geometry type:  POINT
dimension:      XY
bbox:           xmin: 1 ymin: 1 xmax: 3 ymax: 3
CRS:            NA
    a.x a.y        geom
1     1  NA POINT (1 1)
2     2  12 POINT (2 2)
2.1   2  13 POINT (2 2)
3     3  14 POINT (3 3)
> st_join(a, b, left = FALSE)
Simple feature collection with 3 features and 2 fields
geometry type:  POINT
dimension:      XY
bbox:           xmin: 2 ymin: 2 xmax: 3 ymax: 3
CRS:            NA
    a.x a.y        geom
2     2  12 POINT (2 2)
2.1   2  13 POINT (2 2)
3     3  14 POINT (3 3)
> # st_join, largest = TRUE:
> nc <- st_transform(st_read(system.file("shape/nc.shp", package="sf")), 2264)
Reading layer `nc' from data source `/home/edzer/git/sf.Rcheck/sf/shape/nc.shp' using driver `ESRI Shapefile'
Simple feature collection with 100 features and 14 fields
geometry type:  MULTIPOLYGON
dimension:      XY
bbox:           xmin: -84.32385 ymin: 33.88199 xmax: -75.45698 ymax: 36.58965
geographic CRS: NAD27
> gr = st_sf(
+     label = apply(expand.grid(1:10, LETTERS[10:1])[,2:1], 1, paste0, collapse = " "),
+ 	geom = st_make_grid(st_as_sfc(st_bbox(nc))))
> gr$col = sf.colors(10, categorical = TRUE, alpha = .3)
> # cut, to check, NA's work out:
> gr = gr[-(1:30),]
> st_join(nc, gr, largest = TRUE)
Simple feature collection with 100 features and 16 fields
geometry type:  MULTIPOLYGON
dimension:      XY
bbox:           xmin: 406262.2 ymin: 48374.87 xmax: 3052887 ymax: 1044158
projected CRS:  NAD83 / North Carolina (ftUS)
First 10 features:
    AREA PERIMETER CNTY_ CNTY_ID        NAME  FIPS FIPSNO CRESS_ID BIR74 SID74
1  0.114     1.442  1825    1825        Ashe 37009  37009        5  1091     1
2  0.061     1.231  1827    1827   Alleghany 37005  37005        3   487     0
3  0.143     1.630  1828    1828       Surry 37171  37171       86  3188     5
4  0.070     2.968  1831    1831   Currituck 37053  37053       27   508     1
5  0.153     2.206  1832    1832 Northampton 37131  37131       66  1421     9
6  0.097     1.670  1833    1833    Hertford 37091  37091       46  1452     7
7  0.062     1.547  1834    1834      Camden 37029  37029       15   286     0
8  0.091     1.284  1835    1835       Gates 37073  37073       37   420     0
9  0.118     1.421  1836    1836      Warren 37185  37185       93   968     4
10 0.124     1.428  1837    1837      Stokes 37169  37169       85  1612     1
   NWBIR74 BIR79 SID79 NWBIR79 label       col                       geometry
1       10  1364     0      19  A  4 #fb80724d MULTIPOLYGON (((1270814 913...
2       10   542     3      12  A  4 #fb80724d MULTIPOLYGON (((1340555 959...
3      208  3616     6     260  A  5 #80b1d34d MULTIPOLYGON (((1570591 910...
4      123   830     2     145  A 10 #bc80bd4d MULTIPOLYGON (((2881208 948...
5     1066  1606     3    1197  A  8 #fccde54d MULTIPOLYGON (((2525701 911...
6      954  1838     5    1237  A  9 #d9d9d94d MULTIPOLYGON (((2665113 911...
7      115   350     2     139  A 10 #bc80bd4d MULTIPOLYGON (((2881208 948...
8      254   594     2     371  A  9 #d9d9d94d MULTIPOLYGON (((2717989 951...
9      748  1190     2     844  A  8 #fccde54d MULTIPOLYGON (((2203890 914...
10     160  2038     5     176  A  5 #80b1d34d MULTIPOLYGON (((1697624 911...
Warning message:
attribute variables are assumed to be spatially constant throughout all geometries 
> 
> # rbind:
> x = st_sf(a = 1:2, geom = st_sfc(list(st_point(0:1), st_point(0:1)), crs = 4326))
> rbind(x, x, x)
Simple feature collection with 6 features and 1 field
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 1 xmax: 0 ymax: 1
geographic CRS: WGS 84
  a        geom
1 1 POINT (0 1)
2 2 POINT (0 1)
3 1 POINT (0 1)
4 2 POINT (0 1)
5 1 POINT (0 1)
6 2 POINT (0 1)
> nc2 = rbind(nc[1:50, ], nc[51:100, ])
> all.equal(nc, nc2)
[1] TRUE
> 
> # st_sample:
> suppressWarnings(RNGversion("3.5.3"))
> set.seed(131)
> options(digits=6)
> x = st_sfc(st_polygon(list(rbind(c(0,1),c(90,1),c(90,90),c(0,90),c(0,1)))), crs = st_crs(4326))
> #if (sf_extSoftVersion()["proj.4"] >= "4.9.0")
>   (p <- st_sample(x, 10))
although coordinates are longitude/latitude, st_relate_pattern assumes that they are planar
although coordinates are longitude/latitude, st_intersects assumes that they are planar
although coordinates are longitude/latitude, st_intersects assumes that they are planar
Geometry set for 10 features 
geometry type:  POINT
dimension:      XY
bbox:           xmin: 11.2448 ymin: 3.16385 xmax: 82.3451 ymax: 60.1703
geographic CRS: WGS 84
First 5 geometries:
POINT (18.5793 25.2416)
POINT (11.2448 20.8596)
POINT (26.3946 60.1703)
POINT (33.8202 19.146)
POINT (76.1712 32.1029)
>   p <- st_sample(x[[1]], 10) # sfg method
> x = st_sfc(st_polygon(list(rbind(c(0,0),c(90,0),c(90,90),c(0,90),c(0,0))))) # NOT long/lat:
> p <- st_sample(x, 10)
> x = st_sfc(st_polygon(list(rbind(c(-180,-90),c(180,-90),c(180,90),c(-180,90),c(-180,-90)))),
+  crs=st_crs(4326))
> if (sf_extSoftVersion()["proj.4"] >= "4.9.0") # lwgeom breaks on this
+   (p <- st_sample(x, 10))
although coordinates are longitude/latitude, st_relate_pattern assumes that they are planar
although coordinates are longitude/latitude, st_intersects assumes that they are planar
although coordinates are longitude/latitude, st_intersects assumes that they are planar
Geometry set for 10 features 
geometry type:  POINT
dimension:      XY
bbox:           xmin: -148.196 ymin: -51.7878 xmax: 164.644 ymax: 79.4122
geographic CRS: WGS 84
First 5 geometries:
POINT (-148.196 72.3139)
POINT (143.075 54.9954)
POINT (-143.112 25.5397)
POINT (-89.2235 -18.4014)
POINT (-104.569 79.4122)
> pt = st_multipoint(matrix(1:20,,2))
> st_sample(p, 3)
Geometry set for 1 feature 
geometry type:  MULTIPOINT
dimension:      XY
bbox:           xmin: -112.989 ymin: -51.7878 xmax: 53.2886 ymax: 79.4122
CRS:            NA
MULTIPOINT ((-104.569 79.4122), (-112.989 30.33...
> try(st_sample(p, 3.3))
Error in st_sample.sfc(p, 3.3) : size should be an integer
> ls = st_sfc(st_linestring(rbind(c(0,0),c(0,1))),
+  st_linestring(rbind(c(0,0),c(.1,0))),
+  st_linestring(rbind(c(0,1),c(.1,1))),
+  st_linestring(rbind(c(2,2),c(2,2.00001))))
> st_sample(ls, 80)
Geometry set for 4 features  (with 1 geometry empty)
geometry type:  MULTIPOINT
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 0.0875734 ymax: 1
CRS:            NA
MULTIPOINT ((0 0.40852), (0 0.415459), (0 0.999...
MULTIPOINT ((0.00320272 0), (0.0642832 0), (0.0...
MULTIPOINT ((0.0147658 1), (0.0201495 1), (0.04...
MULTIPOINT EMPTY
> st_sample(nc[1:2,], size = c(10,20))
Geometry set for 30 features 
geometry type:  POINT
dimension:      XY
bbox:           xmin: 1200740 ymin: 942550 xmax: 1425790 ymax: 1036860
projected CRS:  NAD83 / North Carolina (ftUS)
First 5 geometries:
POINT (1228809 1011255)
POINT (1233859 956721)
POINT (1260039 1007661)
POINT (1271553 1021505)
POINT (1222552 980314)
> # try with LINES, LongLat, should generate a warning:
> nc[1:2,] %>% st_transform(4326) %>% st_cast("MULTILINESTRING") %>% st_sample(size = c(10,20))
although coordinates are longitude/latitude, st_sample assumes that they are planar
although coordinates are longitude/latitude, st_sample assumes that they are planar
Geometry set for 2 features 
geometry type:  MULTIPOINT
dimension:      XY
bbox:           xmin: -81.7408 ymin: 36.2735 xmax: -80.9188 ymax: 36.5821
geographic CRS: WGS 84
MULTIPOINT ((-81.3308 36.5204), (-81.5205 36.58...
MULTIPOINT ((-81.3127 36.4804), (-81.2655 36.43...
> st_sample(ls, 80, type = "regular")
Geometry set for 4 features  (with 1 geometry empty)
geometry type:  MULTIPOINT
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 0.0987443 ymax: 1
CRS:            NA
MULTIPOINT ((0 0.00373519), (0 0.0187353), (0 0...
MULTIPOINT ((0.00874356 0), (0.0237437 0), (0.0...
MULTIPOINT ((0.0137444 1), (0.0287446 1), (0.04...
MULTIPOINT EMPTY
> p_sample = lapply(1:10, function(i) st_sample(nc[i, ], 100, exact = FALSE))
> lengths(p_sample)
 [1] 106 107 111  92  95 122  86 105 105 101
> p_sample_exact = lapply(1:10, function(i) st_sample(nc[i, ], 100, exact = TRUE))
> lengths(p_sample_exact)
 [1] 100 100 100 100 100 100 100 100 100 100
> #plot(nc$geometry[1])
> #plot(p_sample[[1]], add = TRUE)
> #plot(p_sample_exact[[1]], add = TRUE)
> 
> #class(st_bind_cols(nc, as.data.frame(nc)[1:3]))
> class(dplyr::bind_cols(nc, as.data.frame(nc)[1:3]))
[1] "sf"         "data.frame"
> class(rbind(nc, nc))
[1] "sf"         "data.frame"
> class(cbind(nc, nc))
[1] "sf"         "data.frame"
> 
> x = st_sfc(st_point(0:1), st_point(2:3))
> x[c(NA,1,NA,2,NA)]
Geometry set for 5 features  (with 3 geometries empty)
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 1 xmax: 2 ymax: 3
CRS:            NA
POINT EMPTY
POINT (0 1)
POINT EMPTY
POINT (2 3)
POINT EMPTY
> 
> # jitter
> pts = st_centroid(st_geometry(nc))
> plot(pts)
> plot(st_jitter(pts, .05), add = TRUE, col = 'red')
> plot(st_geometry(nc))
> plot(st_jitter(st_geometry(nc), factor = .01), add = TRUE, col = '#ff8888')
> st_jitter(st_sfc(st_point(0:1)), amount = .1)
Geometry set for 1 feature 
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0.0102133 ymin: 1.06372 xmax: 0.0102133 ymax: 1.06372
CRS:            NA
POINT (0.0102133 1.06372)
> 
> # st_bbox:
> library(sp)
> demo(meuse, ask = FALSE, echo = FALSE)
> suppressWarnings(st_bbox(meuse))
  xmin   ymin   xmax   ymax 
178605 329714 181390 333611 
> crs = suppressWarnings(st_crs(meuse))
> library(raster)

Attaching package: 'raster'

The following object is masked from 'package:dplyr':

    select

> suppressWarnings(st_bbox(raster(meuse.grid)))
  xmin   ymin   xmax   ymax 
178440 329600 181560 333760 
> st_bbox(extent(raster()))
xmin ymin xmax ymax 
-180  -90  180   90 
> 
> # st_to_s2
> if (FALSE) { # stops working with GDAL 2.3.0 / PROJ 5.0.1:
+  x = sf:::st_to_s2(nc)
+  x1 = st_geometry(x)
+  cc = st_coordinates(x1)
+  summary(sqrt(cc[,1]^2+cc[,2]^2+cc[,3]^2))
+ }
> 
> # check_ring_dir
> m = rbind(c(0,0), c(0,1), c(1,1), c(1,0), c(0,0))
> mi = m[nrow(m):1,]
> pol = st_polygon(list(m * 10, m + .5, mi + 1.5, mi + 3.5, m + 5, mi + 6.5))
> st_sfc(pol)
Geometry set for 1 feature 
geometry type:  POLYGON
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 10 ymax: 10
CRS:            NA
POLYGON ((0 0, 0 10, 10 10, 10 0, 0 0), (0.5 0....
> x = st_sfc(pol, check_ring_dir=TRUE)
> y = st_sf(a = 1, geom = st_sfc(pol), check_ring_dir=TRUE)
> str(x)
sfc_POLYGON of length 1; first list element: List of 6
 $ : num [1:5, 1:2] 0 10 10 0 0 0 0 10 10 0
 $ : num [1:5, 1:2] 0.5 0.5 1.5 1.5 0.5 0.5 1.5 1.5 0.5 0.5
 $ : num [1:5, 1:2] 1.5 1.5 2.5 2.5 1.5 1.5 2.5 2.5 1.5 1.5
 $ : num [1:5, 1:2] 3.5 3.5 4.5 4.5 3.5 3.5 4.5 4.5 3.5 3.5
 $ : num [1:5, 1:2] 5 5 6 6 5 5 6 6 5 5
 $ : num [1:5, 1:2] 6.5 6.5 7.5 7.5 6.5 6.5 7.5 7.5 6.5 6.5
 - attr(*, "class")= chr [1:3] "XY" "POLYGON" "sfg"
> x = st_sfc(st_polygon(), st_polygon(), check_ring_dir=TRUE)
> str(x)
sfc_POLYGON of length 2; first list element:  list()
 - attr(*, "class")= chr [1:3] "XY" "POLYGON" "sfg"
> # empty ring/zero area:
> x = st_sfc(st_polygon(list(m[c(1,3,1),])), check_ring_dir=TRUE)
> 
> mp = st_multipolygon(list(pol, pol))
> try(x <- st_sfc(mp, st_polygon(), check_ring_dir=TRUE))
Error in check_ring_dir(lst) : 
  check_ring_dir: not supported for class sfc_GEOMETRY
> x <- st_sfc(mp, pol) %>% st_cast("MULTIPOLYGON") %>% st_sfc(check_ring_dir=TRUE)
> x
Geometry set for 2 features 
geometry type:  MULTIPOLYGON
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 10 ymax: 10
CRS:            NA
MULTIPOLYGON (((0 0, 10 0, 10 10, 0 10, 0 0), (...
MULTIPOLYGON (((0 0, 10 0, 10 10, 0 10, 0 0), (...
> str(x)
sfc_MULTIPOLYGON of length 2; first list element: List of 2
 $ :List of 6
  ..$ : num [1:5, 1:2] 0 10 10 0 0 0 0 10 10 0
  ..$ : num [1:5, 1:2] 0.5 0.5 1.5 1.5 0.5 0.5 1.5 1.5 0.5 0.5
  ..$ : num [1:5, 1:2] 1.5 1.5 2.5 2.5 1.5 1.5 2.5 2.5 1.5 1.5
  ..$ : num [1:5, 1:2] 3.5 3.5 4.5 4.5 3.5 3.5 4.5 4.5 3.5 3.5
  ..$ : num [1:5, 1:2] 5 5 6 6 5 5 6 6 5 5
  ..$ : num [1:5, 1:2] 6.5 6.5 7.5 7.5 6.5 6.5 7.5 7.5 6.5 6.5
  ..- attr(*, "class")= chr [1:3] "XY" "POLYGON" "sfg"
 $ :List of 6
  ..$ : num [1:5, 1:2] 0 10 10 0 0 0 0 10 10 0
  ..$ : num [1:5, 1:2] 0.5 0.5 1.5 1.5 0.5 0.5 1.5 1.5 0.5 0.5
  ..$ : num [1:5, 1:2] 1.5 1.5 2.5 2.5 1.5 1.5 2.5 2.5 1.5 1.5
  ..$ : num [1:5, 1:2] 3.5 3.5 4.5 4.5 3.5 3.5 4.5 4.5 3.5 3.5
  ..$ : num [1:5, 1:2] 5 5 6 6 5 5 6 6 5 5
  ..$ : num [1:5, 1:2] 6.5 6.5 7.5 7.5 6.5 6.5 7.5 7.5 6.5 6.5
  ..- attr(*, "class")= chr [1:3] "XY" "POLYGON" "sfg"
 - attr(*, "class")= chr [1:3] "XY" "MULTIPOLYGON" "sfg"
> 
> x = st_sfc(st_linestring(rbind(c(-179,0),c(179,0))), crs = 4326)
> st_wrap_dateline(st_sf(a = 1, geometry = x))
Simple feature collection with 1 feature and 1 field
geometry type:  MULTILINESTRING
dimension:      XY
bbox:           xmin: -180 ymin: 0 xmax: 180 ymax: 0
geographic CRS: WGS 84
  a                       geometry
1 1 MULTILINESTRING ((-179 0, -...
> st_wrap_dateline(x)
Geometry set for 1 feature 
geometry type:  MULTILINESTRING
dimension:      XY
bbox:           xmin: -180 ymin: 0 xmax: 180 ymax: 0
geographic CRS: WGS 84
MULTILINESTRING ((-179 0, -180 0), (180 0, 179 0))
> st_wrap_dateline(x[[1]])
MULTILINESTRING ((-179 0, -180 0), (180 0, 179 0))
> 
> geo <- c("{\"geodesic\":true,\"type\":\"Point\",\"coordinates\":[-118.68152563269095,36.43764870908927]}",
+          "{\"geodesic\":true,\"type\":\"Point\",\"coordinates\":[-118.67408758213843,36.43366018922779]}",
+          "{\"geodesic\":true,\"type\":\"Point\",\"coordinates\":[-118.67708346361097,36.44208638659282]}",
+          "{\"geodesic\":true,\"type\":\"Point\",\"coordinates\":[-118.67886661944996,36.44110273135671]}",
+          "{\"geodesic\":true,\"type\":\"Point\",\"coordinates\":[-118.68089232041565,36.44173155205561]}")
> st_as_sfc(geo, GeoJSON = TRUE)
Geometry set for 5 features 
geometry type:  POINT
dimension:      XY
bbox:           xmin: -118.682 ymin: 36.4337 xmax: -118.674 ymax: 36.4421
geographic CRS: WGS 84
POINT (-118.682 36.4376)
POINT (-118.674 36.4337)
POINT (-118.677 36.4421)
POINT (-118.679 36.4411)
POINT (-118.681 36.4417)
> st_as_sfc(geo, GeoJSON = TRUE, crs = 4326)
Geometry set for 5 features 
geometry type:  POINT
dimension:      XY
bbox:           xmin: -118.682 ymin: 36.4337 xmax: -118.674 ymax: 36.4421
geographic CRS: WGS 84
POINT (-118.682 36.4376)
POINT (-118.674 36.4337)
POINT (-118.677 36.4421)
POINT (-118.679 36.4411)
POINT (-118.681 36.4417)
> 
> st_as_sfc(st_as_binary(st_sfc(st_point(0:1)))[[1]], crs = 4326)
Geometry set for 1 feature 
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 1 xmax: 0 ymax: 1
geographic CRS: WGS 84
POINT (0 1)
> 
> x = nc
> x$geom = NULL
> class(x)
[1] "sf"         "data.frame"
> 
> st_as_sfc(list(st_point(0:1)), crs = 4326)
Geometry set for 1 feature 
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 1 xmax: 0 ymax: 1
geographic CRS: WGS 84
POINT (0 1)
> 
> # crop:
> box = c(xmin = 0, ymin = 0, xmax = 1, ymax = 1)
> 
> pol = st_sfc(st_buffer(st_point(c(.5, .5)), .65))
> pol_sf = st_sf(a=1, geom=pol)
> 
> st_crop(pol, box)
Geometry set for 1 feature 
geometry type:  POLYGON
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
POLYGON ((1 0.0849179, 0.983044 0.0650651, 0.95...
> st_crop(pol, st_bbox(box))
Geometry set for 1 feature 
geometry type:  POLYGON
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
POLYGON ((1 0.0849179, 0.983044 0.0650651, 0.95...
> st_crop(pol_sf, box)
Simple feature collection with 1 feature and 1 field
geometry type:  POLYGON
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
  a                           geom
1 1 POLYGON ((1 0.0849179, 0.98...
Warning message:
attribute variables are assumed to be spatially constant throughout all geometries 
> st_crop(pol_sf, st_bbox(box))
Simple feature collection with 1 feature and 1 field
geometry type:  POLYGON
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 1 ymax: 1
CRS:            NA
  a                           geom
1 1 POLYGON ((1 0.0849179, 0.98...
Warning message:
attribute variables are assumed to be spatially constant throughout all geometries 
> 
> # new sample methods:
> x = st_sfc(st_polygon(list(rbind(c(0,0),c(90,0),c(90,90),c(0,90),c(0,0))))) # NOT long/lat:
> p <- st_sample(x, 10, type = "regular")
> p <- st_sample(x, 10, type = "hexagonal")
> 
> all.equal(st_drop_geometry(pol_sf), st_set_geometry(pol_sf, NULL))
[1] TRUE
> 
> # https://github.com/r-spatial/sf/issues/1024
> shape1 <-st_sfc(st_polygon(list(rbind(c(0,0), c(1,0), c(3,2), c(2,4), c(1,4), c(0,0)))))
> shape2 <- st_sfc(st_polygon())
> shape3 <- st_sfc(st_polygon())
> 
> shape4 = st_intersection(shape2, shape3) # has zero features
> 
> st_difference(shape1, shape4)
Geometry set for 1 feature 
geometry type:  POLYGON
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 3 ymax: 4
CRS:            NA
POLYGON ((0 0, 1 0, 3 2, 2 4, 1 4, 0 0))
> st_difference(shape4, shape1)
Geometry set for 0 features 
bbox:           xmin: NA ymin: NA xmax: NA ymax: NA
CRS:            NA
> st_sym_difference(shape1, shape4)
Geometry set for 1 feature 
geometry type:  POLYGON
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 3 ymax: 4
CRS:            NA
POLYGON ((0 0, 1 0, 3 2, 2 4, 1 4, 0 0))
> st_union(shape1, shape4)
Geometry set for 1 feature 
geometry type:  POLYGON
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 3 ymax: 4
CRS:            NA
POLYGON ((0 0, 1 0, 3 2, 2 4, 1 4, 0 0))
> st_union(shape4, shape1)
Geometry set for 1 feature 
geometry type:  POLYGON
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 3 ymax: 4
CRS:            NA
POLYGON ((0 0, 1 0, 3 2, 2 4, 1 4, 0 0))
> 
> # transform empty:
> tr = st_sf(geom=st_sfc()) %>% st_set_crs(3587) %>% st_transform(4326)
> 
> # NA values are converted to empty; #1114:
> x <- data.frame(name=LETTERS)
> y <- data.frame(name=LETTERS[1:13], letters[14:26])
> y$geometry <- st_sfc(st_point(c(0,0)))
> y <- st_sf(y)
> out = merge(x, y, all.x=TRUE)
> class(out)
[1] "data.frame"
> 
> st_as_sf(st_sfc(st_point(0:1)))
Simple feature collection with 1 feature and 0 fields
geometry type:  POINT
dimension:      XY
bbox:           xmin: 0 ymin: 1 xmax: 0 ymax: 1
CRS:            NA
            x
1 POINT (0 1)
> 
> proc.time()
   user  system elapsed 
  6.032   0.160   6.217