r-spatial / sf Goto Github PK

View Code? Open in Web Editor NEW

1.3K 64.0 291.0 131.05 MB

Simple Features for R

Home Page: https://r-spatial.github.io/sf/

License: Other

R 64.15% C++ 30.21% Shell 0.01% C 0.24% Dockerfile 3.57% M4 1.79% Awk 0.03%

r spatial gdal geos proj r-package rstats

sf's Introduction

Simple Features for R

A package that provides simple features access for R.

Blogs, links • Cheatsheet • Installing • Contributing • Acknowledgment • How to cite

Package sf:

represents simple features as records in a data.frame or tibble with a geometry list-column
represents natively in R all 17 simple feature types for all dimensions (XY, XYZ, XYM, XYZM)
interfaces to GEOS for geometrical operations on projected coordinates, and (through R package s2) to s2geometry for geometrical operations on ellipsoidal coordinates
interfaces to GDAL, supporting all driver options, Date and POSIXct and list-columns
interfaces to PRØJ for coordinate reference system conversion and transformation
uses well-known-binary serialisations written in C++/Rcpp for fast I/O with GDAL and GEOS
reads from and writes to spatial databases such as PostGIS using DBI
is extended by
- lwgeom for selected liblwgeom/PostGIS functions
- stars for raster data, and raster or vector data cubes (spatial time series)
- sfnetworks for geospatial network data

Books, journal articles, blogs, presentations, vignettes, sp-sf wiki

an open access R Journal article summarizes the package
two books: Spatial Data Science: with applications in R, Geocomputation with R
package vignettes: first, second, third, fourth, fifth, sixth, seventh
blog posts: first, second, third, fourth
the original R Consortium ISC proposal, the R Consortium blog post
presentations: rstudio::conf 2018 (video), UseR! 2016
wiki page describing sp-sf migration

Cheatsheet

CC 4.0 BY Ryan Garnett

Installing

Install either from CRAN with:

install.packages("sf")

This will install binary packages on Windows and MacOS, unless you configured R such that it tries to install source packages; in that case, see below.

Install development versions from GitHub with:

library(remotes)
install_github("r-spatial/sf")

Windows

Installing sf from source works under Windows when Rtools is installed.

MacOS

MacOS users are strongly encouraged to install the sf binary packages from CRAN, unless they are familiar with compilers, linking, C++ source code, and homebrew.

The easiest way to install gdal is using Homebrew. Recent versions of Homebrew include a full-featured up-to-date gdal formula, which installs proj and gdal at the same time:

brew install pkg-config
brew install gdal

Once gdal is installed, you will be able to install sf package from source in R. With the current version of proj on homebrew, installation requires additional configuration:

install.packages("sf", type = "source", configure.args = "--with-proj-lib=$(brew --prefix)/lib/")

Or the development version:

library(remotes)
install_github("r-spatial/sf", configure.args = "--with-proj-lib=$(brew --prefix)/lib/")

If you are using sf and rgdal together, it is necessary to install rgdal from source using this configuration:

install.packages("rgdal", type = "source", configure.args = c("--with-proj-lib=$(brew --prefix)/lib/", "--with-proj-include=$(brew --prefix)/include/"))

Alternatively, these instructions explain how to install gdal using kyngchaos frameworks.

For Mac OS 11 Big Sur source install instruction, see here

Linux

For Unix-alikes, GDAL (>= 2.0.1), GEOS (>= 3.4.0) and Proj.4 (>= 4.8.0) are required.

Ubuntu

Dependencies for recent versions of Ubuntu (18.04 and later) are available in the official repositories; you can install them with:

apt-get -y update && apt-get install -y  libudunits2-dev libgdal-dev libgeos-dev libproj-dev

However, to get more up-to-date versions of dependencies such as GDAL, we recommend adding the ubuntugis-unstable PPA to the package repositories and installing them as follows:

sudo add-apt-repository ppa:ubuntugis/ubuntugis-unstable
sudo apt-get update
sudo apt-get install libudunits2-dev libgdal-dev libgeos-dev libproj-dev libsqlite0-dev

Adding this PPA is required for installing sf on older versions of Ubuntu (e.g. Xenial).

Another option, for advanced users, is to install dependencies from source; see e.g. an older Travis config file for hints.

Fedora

The following command installs all required dependencies:

sudo dnf install gdal-devel proj-devel geos-devel sqlite-devel udunits2-devel

Arch

Get gdal, proj, geos and podofo from the main repos, and udunits from the AUR:

pacman -S gdal proj geos arrow podofo
yay/pacaur/yaourt/whatever -S udunits

`renv` or `conda`

There are several reports that sf fails to install as a source package when R is used with renv, or when R is installed in a conda environment. If you experience this, please do not raise an issue here, but

try to sort this out with the renv developers or the conda maintainers, or
try to use binary installs of the sf package, e.g. from r2u, or the Posit package manager

Other

To install on Debian, the rocker geospatial Dockerfiles may be helpful. Ubuntu Dockerfiles are found here.

Multiple GDAL, GEOS and/or PROJ versions on your system

If you use dynamic linking (installation from source) and have multiple versions of these libraries installed (e.g. one from ubuntugis-unstable, another installed from source in /usr/local/lib) then this will in general not work, even when setting LD_LIBRARY_PATH manually. See here for the reason why.

lwgeom

Functions and methods that require liblwgeom, including ellipsoidal (not spherical or Eucledian) metrics (area, distances), are provide by and used from lwgeom, which is also on CRAN.

Contributing

Contributions of all sorts are most welcome, issues and pull requests are the preferred ways of sharing them.
When contributing pull requests, please adhere to the package style (in package code use = rather than <-; don't change indentation; tab stops of 4 spaces are preferred)
This project is released with a Contributor Code of Conduct. By participating in this project you agree to abide by its terms.

How to cite

Package sf can be cited as:

Edzer Pebesma, 2018. Simple Features for R: Standardized Support for Spatial Vector Data. The R Journal 10:1, 439-446.
Pebesma, E.; Bivand, R. (2023). Spatial Data Science: With Applications in R (1st ed.). 314 pages. Chapman and Hall/CRC.

Acknowledgment

This project gratefully acknowledges financial support from the

sf's People

Contributors

Stargazers

Watchers

Forkers

jefferis kendonb bleutner etiennebr jimhester nuest greenstat epsimatic88 cpsievert dheerajpatta loicdtx thk686 barryrowlingson bengraeler astroidex demorenoc 13671359855 odeleongt hmenag1 bbest ecor danovando pythianconsultants lbusett zhaoxiaohe yutannihilation ateucher nemochina2008 tim-salabim sbakuleva gergness zedseayou geocfllc hughjonesd matthieustigler wangrandk buyuhuang grseb9s fototo dshkol theresa-gao rsbivand johnbaums cmcaine valexandre christophertull robinlovelace nikolayvoronchikhin bhavmakwana marek1914 choisy 3014904018 dbaston statnmap anhnguyendepocen tghoward djrondon uribo jbollow jahonamir sherinsh humbertosubiza khondula ejrivera6 jmlondon radovankavicky gapdata ekarsten dpprdan soonho2 guzhongru haozhu233 tylerlittlefield ryfx monsieurxu gisma 10square sanderdevisscher linzhu2017 gravitytrope ysosirius rnuske angela-li lambdamoses vfonov thomasp85 jmarshallnz hhy5277 robitalec nikolai-b arielhirschhorn luigibernardi hannahbus vlulla adrianauscanga stewartmacdonald doudatadog symbolixau mikejohnson51 chihyulai

sf's Issues

Missing sf Windows binary proj-datumgrid files

This example taken from the PROJ documentation fails:

p1 = st_point(c(-111, 50))
sfc = st_sfc(p1,
 crs = "+proj=longlat +ellps=clrk66 +nadgrids=ntv1_can.dat")
st_transform(sfc, "+proj=latlong +ellps=GRS80 +datum=NAD83")

on Windows with

> library(sf)
Linking to GEOS 3.5.0, GDAL 2.1.1
> p1 = st_point(c(-111, 50))
> sfc = st_sfc(p1,
+  crs = "+proj=longlat +ellps=clrk66 +nadgrids=ntv1_can.dat")
> st_transform(sfc, "+proj=latlong +ellps=GRS80 +datum=NAD83")
Error code: 6
Error in eval(substitute(expr), envir, enclos) : OGR error

because the ntv1_can.dat file is not in sf/proj. Commonly used files are in proj-datumgrid*.zip, and should be available.

The call goes through ogr/ogrct.cpp where attempts are made to handle error messages, obviously not very effectively on the OGR/GDAL side.

When the datum grid is present, we see:

> library(sf)
Linking to GEOS 3.6.0, GDAL 2.1.2
> p1 = st_point(c(-111, 50))
> sfc = st_sfc(p1,
+  crs = "+proj=longlat +ellps=clrk66 +nadgrids=ntv1_can.dat")
> st_transform(sfc, "+proj=latlong +ellps=GRS80 +datum=NAD83")
Geometry set for 1 feature 
geometry type:  POINT
dimension:      XY
bbox:           xmin: -111.0008 ymin: 50.00003 xmax: -111.0008 ymax: 50.00003
epsg (SRID):    NA
proj4string:    +proj=latlong +ellps=GRS80 +datum=NAD83
POINT(-111.000820033329 50.0000308342618) 
> c(-(111+(2.952/3600)), 50+(0.111/3600))
[1] -111.00082   50.00003

st_read and requiring layer name

This is a small pet peeve of mine with regard to rgdal - is there a compelling reason to continue to force the user to provide the layer name when reading in data?

There isn't a technical reason within gdal to require this and in my experience most geometry files I've worked with only ever have one layer. It seems to me a more sensible approach is to make the layer name an optional argument and in the case of a single layer to just automatically return that object.

If more than one layer are present then the current error could be thrown, or alternatively throw a warning and return a list containing all layers.

Int64, package bit64

GDAL 2.x can read Int64 values, but (base) R can't store them. Right now, sf converts these to double, and warns when information gets lost (values are over 2^53). Another approach would be to return them as integer64, an S3 type provided by package bit64.

What about the M component?

There is some discussion on gdal-dev about the M-component (measure) when points are not just 2D. There is an open issue on the handling of M as an attribute in non-POINT features. Z and M are present in OGC.

unknown wkbType dim in switch

Why is st_read_db dependent on the PostgreSQL driver used for the connection:

# using <PostgreSQLDriver> from package RPostgreSQL
conn <- DBI::dbConnect(RPostgreSQL::PostgreSQL(), 
    dbname = "postgis")
st_read_db(conn, query = "select geom from plots limit 10;")

Simple feature collection with 10 features and 0 fields
geometry type:  POINT
dimension:      XY
bbox:           xmin: 1802914 ymin: 1076334 xmax: 1907419 ymax: 1334777
epsg (SRID):    3175
proj4string:    +proj=aea +lat_1=42.122774 +lat_2=49.01518 +lat_0=45.568977 +lon_0=-83.248627 +x_0=1000000 +y_0=1000000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
precision:      double (default; no precision model)
                             geom
1  POINT(1849836.61661095 1076...
2  POINT(1849806.67270222 1076...
3  POINT(1875028.97521343 1116...
4  POINT(1874767.55838702 1116...
5  POINT(1803250.35520253 1334...
6  POINT(1802913.8109178 13347...
7  POINT(1907418.52662745 1223...
8  POINT(1840845.46786735 1220...
9  POINT(1836423.52704828 1226...
10 POINT(1836041.74915185 1226...

# using <PqDriver> from package RPostgres
conn <- DBI::dbConnect(RPostgres::Postgres(), 
    dbname = "postgis")
st_read_db(conn, query = "select geom from plots limit 10;")

Error in eval(substitute(expr), envir, enclos) : 
  unknown wkbType dim in switch

However, using this driver works if I use the results from a previous query:

library(data.table)
conn <- DBI::dbConnect(RPostgres::Postgres(), 
    dbname = "postgis")
plots <- data.table(DBI::dbGetQuery(conn,
    dplyr::build_sql("
        SELECT   id_pe, geom
        FROM     plots
        LIMIT    10;")))
plots[, sfc := st_as_sfc(structure(geom, class = "WKB"), EWKB = TRUE, 
    crs = "+init:epsg=3175")]
plots[, str(sfc)]

Is something missing in the code for st_read_db (or absent from the doc) or does something differ in the way each driver return geometries (or I miss extra options when calling st_read_db with driver PqDriver)?

Building OSX GDAL/PROJ.4/GEOS and ext-deps for GDAL drivers

I guess we still need --enable-shared=no in PROJ, and in GDAL --with-static-proj4= and --enable-shared=no. I believe CRAN had a good view of this, but I don't know where their logs are. GDAL drivers are important, and should also be built static (I guess).

Set CRS in a pipe - idea

Below you can find my broad idea how it should work. Basically, I think there should be a function to set a crs in a pipe, such as st_crs(sf_object, crs).

library('sf')
library('dplyr')
filepath <- system.file("shape/nc.shp", package="sf")
nc <- st_read(filepath) %>% st_crs(., "4326")
nc <- st_read(filepath) %>% st_crs(., "+init=epsg:4326")

st_transform will not handle +proj=obtrans

Just a placeholder - see the obtrans notes in rgdal::project and rgdal::spTransform. For its own reasons, proj swaps from and to crs here.

Please consider implementing st_makevalid

I've searched many times over the years for an R version of postgis' st_makevalid() function (via liblwgeom) and wonder if this would be a good time and place for it. Thanks for all your hard work!

Slowdown when printing a simple feature

This script works fine when loading the data, but seems to stall when you want to print it:

download.file("http://download.geofabrik.de/europe/netherlands-latest-free.shp.zip", destfile = "/tmp/netherlands-latest-free.shp.zip")
unzip("/tmp/netherlands-latest-free.shp.zip", exdir = "/tmp/")
res = sf::st_read(dsn = "/tmp/", "gis.osm_roads_free_1") # 1.2 gb!
res
# head(res) # not tried yet as previous command stalled R

GDAL 2.1.0 out

Please note:

https://lists.osgeo.org/pipermail/gdal-dev/2016-May/044242.html

http://trac.osgeo.org/gdal/wiki/Release/2.1.0-News

https://trac.osgeo.org/gdal/wiki/rfc61_support_for_measured_geometries

https://trac.osgeo.org/gdal/wiki/rfc59.1_utilities_as_a_library

https://trac.osgeo.org/gdal/wiki/rfc48_geographical_networks_support

no definition of "scale model for coordinates"

I'm not sure - is a standard thing and I just don't know what it is?

stringsAsFactors

Would it be best to set this to FALSE in ST_read?

df = data.frame(row.names = ids, apply(do.call(rbind, f), 2, unlist), stringsAsFactors = FALSE)

Arguably this could be user-controlled, but factors is an R feature so it seems best to stick to straight character. Happy to do the PR if you want the user-option control.

Not returning 1st layer or message for osm data

In relation of #39, please test this reproducible example, called test.R:

u = "http://overpass-api.de/api/interpreter?data=%2F*%0AThis%20has%20been%20generated%20by%20the%20overpass-turbo%20wizard.%0AThe%20original%20search%20was%3A%0A%E2%80%9Cleisure%3Dpark%E2%80%9D%0A*%2F%0A%5Bout%3Axml%5D%2F*fixed%20by%20auto%20repair*%2F%5Btimeout%3A25%5D%3B%0A%2F%2F%20gather%20results%0A%28%0A%20%20%2F%2F%20query%20part%20for%3A%20%E2%80%9Cleisure%3Dpark%E2%80%9D%0A%20%20node%5B%22leisure%22%3D%22park%22%5D%2853.794871619441025%2C-1.5554237365722654%2C53.81063582424263%2C-1.5232372283935547%29%3B%0A%20%20way%5B%22leisure%22%3D%22park%22%5D%2853.794871619441025%2C-1.5554237365722654%2C53.81063582424263%2C-1.5232372283935547%29%3B%0A%20%20relation%5B%22leisure%22%3D%22park%22%5D%2853.794871619441025%2C-1.5554237365722654%2C53.81063582424263%2C-1.5232372283935547%29%3B%0A%29%3B%0A%2F%2F%20print%20results%0Aout%20meta%3B%2F*fixed%20by%20auto%20repair*%2F%0A%3E%3B%0Aout%20meta%20qt%3B%2F*fixed%20by%20auto%20repair*%2F"
download.file(url = u, destfile = "test.osm")
psf = st_read("test.osm")
# user: expecting spatial data
psf # user: oh...
psf_p = st_read("test.osm", "points")
psf_o = st_read("test.osm", "multipolygons")
knitr::spin("/tmp/test.R")

I get no message saying it's returning non-spatial data - seems that functionality is not implemented for osm data yet. Incidentally I think in the process of this test I've found an issue with how GDAL reads OSM data - cc @mpadge.

st_write arguments

Similar to #39, but about st_write().
Now, st_read() accepts filename, but to write a file using st_write there is a need to specify layer argument:

st_write(nc, "nc.shp")
Error in CPL_write_ogr(obj, dsn, layer, driver, as.character(dataset_options),  : 
  argument "layer" is missing, with no default

Is it possible to add this new functionality also for st_write()?

dataset creation options missing

st_read now has an options argument, which accepts layer creation options, but not a way to specify dataset creation options; see #26

saving sf objects with geometry type: GEOMETRY

Is there a recommended file format for storing sf objects that have geometry type: GEOMETRY?

It seems like both .shp and .gpkg return the following error:

GEOMETRY :
Error in eval(substitute(expr), envir, enclos) : unknown type!

A reproducible test:

# Polygon geometry: pol

p1 <- rbind(c(0,0), c(1,0), c(3,2), c(2,4), c(1,4), c(0,0))
p2 <- rbind(c(1,1), c(1,2), c(2,2), c(1,1))
pol <-st_polygon(list(p1,p2))
# 
# Multi-polygon geometry: mpol

p3 <- rbind(c(3,0), c(4,0), c(4,1), c(3,1), c(3,0))
p4 <- rbind(c(3.3,0.3), c(3.8,0.3), c(3.8,0.8), c(3.3,0.8), c(3.3,0.3))[5:1,]
p5 <- rbind(c(3,3), c(4,2), c(4,3), c(3,3))
mpol <- st_multipolygon(list(list(p1,p2), list(p3,p4), list(p5)))

# sf object w/ polygon: sf_pol

sf_pol <- 
        data.frame(VAR1 = 'foo',VAR2 = 'bar') %>% 
        mutate(geom = st_sfc(pol)) %>% 
        st_as_sf()

# sf object w/ multi-polygon: sf_mpol

sf_mpol <-
        data.frame(VAR1 = 'foo',VAR2 = 'bar') %>% 
        mutate(geom = st_sfc(mpol)) %>% 
        st_as_sf()

# combined sf's, geometry type: GEOMETRY

sf_geometry <- rbind(sf_pol,sf_mpol)

sf_geometry

# attempt to write the sf object (currently both fail)
sf_geometry %>% 
        st_write(dsn = '.',layer = 'test',driver = 'ESRI Shapefile')

sf_geometry %>% 
        st_write(dsn = 'test.gpkg',driver = 'GPKG')

[Request] support for raster operations

Does Team sfr plan to also overhaul the raster package? That package has a lot of good features, but there is likely value in having a common interface for all spatial objects that center around data.frames or lists of data.frames.

extract, for example, might go something like this:

Suppose extract(raster_brick, MULTIPOLYGON_list, weights = TRUE) returns a list of data.frames (or tibbles if tibble is loaded) each of which contains a column of extracted values from the RasterBrick (named value), and a column of weights (adding to 1) representing area overlap (named weight).

Then if mp is an sf with a MULTIPOLYGON column named geometry,

mp <- mp %>%
  mutate(raster_band = extract(raster_brick, geometry, weights = TRUE))

would add a column of extracted values as a list of data.frames to mp.

Then, one could do something like this using purrr::map_dbl:

mp <- mp %>%
  mutate(raster_band = extract(raster_brick, geometry, weights = TRUE)) %>%
  mutate(raster_mean = raster_band %>% map_dbl(function(.data){
    sum(.data$weight*.data$value)
  }))

which would add the weighted mean for each row in the original sf object. Alternatively, using tidyr::unnest(), tidyr::nest(), as long as there isn't too much data:

mp <- mp %>%
  mutate(raster_band = extract(raster_brick, geometry, weights = TRUE)) %>%
  unnest(raster_band) %>% 
  group_by(polygon_id) %>% # uniquely identifies polygons
  summarize(raster_mean = sum(weight*value))

would achieve something similar, but would drop some columns.

CRS with added +init=epsg: slot if known

The thread at

https://lists.osgeo.org/pipermail/gdal-dev/2016-May/044474.html

shows that interoperability is harmed when only the p4s is used, not the epsg code (which has authority tags). Could we add an epsg code slot to CRS, default NA, but used when exporting data to resolve to a WKT representation with authority tags?

Another WKT pkg

you mention the wkb package in the readme - i maintain another package for working with WKT data https://github.com/ropensci/wellknown just in case you hadn't seen it and it's of any help

? feather

Do we need to accommodate:

https://blog.rstudio.org/2016/03/29/feather/

Can't connect to PostGIS using st_read()

I can't connect to my PostGIS database using the following command (the same dsn works with ogrinfo):

pol <- st_read("PG:dbname='gis' host='localhost' port=5432 user='mylogin' password='mypassword'", "p1")

I get the following error: Error in eval(substitute(expr), envir, enclos) : Open failed

However the following query is OK with st_read_db():

dbcon <- dbConnect(dbDriver("PostgreSQL"), dbname = "gis", host = "localhost", user = "mylogin", password = "mypassword")
pol <- st_read_db(dbcon, query = "select * from p1")

XYZ shapefiles from worldwind

These look like a good set for testing on:

http://dds.cr.usgs.gov/srtm/version2_1/SWBD/SWBDwest/

found on the Worldwind forum: http://forum.worldwindcentral.com/showthread.php?26821-Sample-3-D-Shapefile

mixed topology files - example GPKG

I was putting together a GPKG with areas, lines, and points but got sidelined by Manifold's export of closed lines - it seems they are written and/or interpreted as polygons, so I have to do some more to make this useable in sf.

Still, it's an example GPKG file.
us_counties_states.zip

library(USAboundaries)
library(sf)

## warns but still works 
# Warning message:
#   In CPL_write_ogr(obj, dsn, layer, driver, as.character(options),  :
#                      GDAL Message 6: dataset us_counties.mif does not support layer creation option FORMAT

st_write(st_as_sf(us_counties()), "us_counties.mif", "us_counties", driver = "MapInfo File", options = "FORMAT=MIF")

## read into Manifold and add lines and points
## EDIT: no Manifold export writes MULTIPOLYGON here 
## dissolve counties by state, add boundaries of the states
## add inner centroids of the states
## points only get state attributes copied
## export to MIF

cs <- st_read("Us_counties_states.mif", "Us_counties_states")
st_write(cs, "us_counties_states.gpkg", "us_counties_states", "GPKG")


cs <- st_read("Us_counties_states.gpkg", "Us_counties_states")
## plotting these the counties-polys and state-points get covered by the state-polys, which is
## not what was intended

SFS is not defined in the document

OGC SFS is not defined as far as I can tell

stringsAsFactors = FALSE not working

Here are a couple of tests to be added in test_read.R

test_that("stringsAsFactors = FALSE produces a data.frame with no factors", {
  nc <-  st_read(system.file("shape/nc.shp", package="sf"), "nc", crs = 4267,
                 stringsAsFactors = FALSE)
  expect_false(any(sapply(nc, class) == "factor"))
})

test_that("stringsAsFactors = TRUE produces a data.frame with factors", {
  nc <-  st_read(system.file("shape/nc.shp", package="sf"), "nc", crs = 4267,
                 stringsAsFactors = TRUE)
  expect_true(any(sapply(nc, class) == "factor"))
})

MapInfo closed linestring interpreted as POLYGON

I've attached two folder data sources "TABBED" and "MIFFED", each contain the same 2 POLYGONs and 2 LINESTRINGs.

linestring_as_polygon.zip

There are four features, in order as:

feature is a POLYGON box
feature is LINESTRING boundary of 1.
feature is a POLYGON box
feature is modfied LINESTRING boundary of 3. with closing vertex moved from origin

ogrinfo sees the same for both folders, namely POLYGON, LINESTRING, POLYGON, LINESTRING

ogrinfo MIFFED Drawing
Had to open data source read-only.
INFO: Open of `MIFFED'
      using driver `MapInfo File' successful.

Layer name: Drawing
Geometry: Unknown (any)
Feature Count: 4
Extent: (-82.000000, -23.500000) - (282.000000, 76.500000)
Layer SRS WKT:
LOCAL_CS["Nonearth",
    UNIT["Meter",1.0]]
ID: Integer (0.0)
OGRFeature(Drawing):1
  ID (Integer) = 1
  Style = BRUSH(fc:#c0c0c0,bc:#808080,id:"mapinfo-brush-2,ogr-brush-0");PEN(w:0pt,c:#808080,id:"mapinfo-pen-2,ogr-pen-0")
  POLYGON ((-82.0 16.5,-82.0 76.5,-22.0 76.5,-22.0 16.5,-82.0 16.5))

OGRFeature(Drawing):2
  ID (Integer) = 2
  Style = PEN(w:1px,c:#000000,id:"mapinfo-pen-2,ogr-pen-0")
  LINESTRING (-82.0 16.5,-82.0 76.5,-22.0 76.5,-22.0 16.5,-82.0 16.5)

OGRFeature(Drawing):3
  ID (Integer) = 3
  Style = BRUSH(fc:#c0c0c0,bc:#808080,id:"mapinfo-brush-2,ogr-brush-0");PEN(w:0pt,c:#808080,id:"mapinfo-pen-2,ogr-pen-0")
  POLYGON ((194.0 -23.5,194.0 62.5,282.0 62.5,282.0 -23.5,194.0 -23.5))

OGRFeature(Drawing):4
  ID (Integer) = 4
  Style = PEN(w:1px,c:#000000,id:"mapinfo-pen-2,ogr-pen-0")
  LINESTRING (189.5 -18.625,194.0 62.5,282.0 62.5,282.0 -23.5,194.0 -23.5)

sf sees 2. as a POLYGON

st_read("TABBED", "Drawing")
## same as below
st_read("MIFFED", "Drawing")
...

precision: double (default; no precision model)
ID geometry
1 1 GEOMETRYCOLLECTION(POLYGON(...
2 2 GEOMETRYCOLLECTION(POLYGON(...
3 3 GEOMETRYCOLLECTION(POLYGON(...
4 4 GEOMETRYCOLLECTION(LINESTRI...

Dplyr compatibility?

I see that working with sf objects with dplyr is listed on the ISC Proposal, but it looks like you haven't gotten around to documenting that yet.

Without going into a vignette-length explanation, could you shed some light on why dplyr functions appear to strip sf objects of their sf class and suggest a way to effectively combine the power of these two tools?

Example:

library(dplyr)
library(sf)

nc <- st_read(system.file("shapes/", package="maptools"), "sids")
## Reading layer sids from data source 
## /Library/Frameworks/R.framework/Versions/3.3/Resources/library/maptools/shapes/ using driver "ESRI Shapefile"
## features:       100
## fields:         14
## converted into: MULTIPOLYGON

class(nc)
## [1] "sf"         "data.frame"

class(nc[nc$AREA > .1,])
## [1] "sf"         "data.frame"

class(dplyr::filter(nc,AREA > .1))
## [1] "data.frame"

Many thanks.

rename sfi to sfg

now, sfi refers to simple feature item. The simple feature however refers to both geometry and attributes, so a better description of a single simple feature's geometry would be sfg, for simple feature geometry.

Calling a collection of these sfgc (c for collection) would be confusing as GEOMETRYCOLLECTION refers to a type of geometry for a single simple feature; I'll stick with sfc, for now (simple feature collection, or column, although not entirely right).

sf build error on CRAN platform macos-x86_64-clang

For CRAN platform r-devel-macos-x86_64-clang I see the following error:

** testing if installed package can be loaded
Error in dyn.load(file, DLLpath = DLLpath, ...) : 
  unable to load shared object '/Users/ripley/R/packages/tests-devel/sf.Rcheck/sf/libs/sf.so':
  dlopen(/Users/ripley/R/packages/tests-devel/sf.Rcheck/sf/libs/sf.so, 6): Symbol not found: _H5T_NATIVE_DOUBLE_g
  Referenced from: /Users/ripley/R/packages/tests-devel/sf.Rcheck/sf/libs/sf.so
  Expected in: flat namespace
 in /Users/ripley/R/packages/tests-devel/sf.Rcheck/sf/libs/sf.so
Error: loading failed
Execution halted
ERROR: loading failed

BDR mentioned that you have not allowed for static builds (on macOS, where it needs to link to hdf5, at least) , but I don't see how to solve this.

Planar index

I just noticed RcppAnnoy, which, because the index is storable after building, might be a candidate for point data.

Extend data.table instead of data.frame?

I've occasionally tried to extend data.frame classes and always given up. I've never found a satisfactory way to store non-trivial things in a data frame column.

See my poorly documented spong package for example: https://github.com/barryrowlingson/spong

BUT... data.table provides a more flexible data grid structure that is very happy to store structured data in columns. Example:

sq=cbind(c(0,1,1,0,0),c(0,0,1,1,0))
polys = list(sq, sq+1, sq+2)
class(polys)="sf"
attr(polys, "CRS")="+init=EPSG:4326"

print.sf=function(x,...){print("Geometry...")}
require(data.table)
d = data.table(polys, name=c("Alpha","Bravo","Charlie"), x=runif(3), y=runif(3))

d now prints like this:

> d
                  polys    name         x         y
1: 0,1,1,0,0,0,0,1,1,0,   Alpha 0.3505690 0.1249515
2: 1,2,2,1,1,1,1,2,2,1,   Bravo 0.5857398 0.5670319
3: 2,3,3,2,2,2,2,3,3,2, Charlie 0.1299579 0.1650072

and the polys column still has its class:

> d$polys
[1] "Geometry..."
`
data.tables are directly compatible with `dplyr` and `ggplot2` (unlike `sp` classes)

d %>% filter(name=="Alpha")
polys name x y
1: 0,1,1,0,0,0,0,1,1,0, Alpha 0.350569 0.1249515
ggplot(d, aes(x=x,y=y)) + geom_point()

Usefully, data table row selection preserves attributes:

> a = d %>% filter(name=="Alpha")
> attr(a$polys,"CRS")
[1] "+init=EPSG:4326"

which is something even R's default subsetting doesn't do - it drops as much as it can including the class of the object:

> class(polys)
[1] "sf"
> class(polys[1])
[1] "list"

and hence you waste your life writing [.sf methods that do little more than restore the attributes that R took away in the first place (see for example getAnywhere("[.POSIXct"))

The only basic thing I can't figure out at the moment is how to identify a geometry column within a data.table. We could add an attribute to a data table, but that gets lost on selection - perhaps the data.table authors might like to help with this:

> attr(d,"geom")="polys"
> attr(d[1,],"geom")
NULL

or we define a superclass of data.table and write some methods for that.

SpatialPolygonsDataTable anyone?

Copying proj-datumgrids to Windows and OSX binary packages

The current proj-datumgrids file should be unpacked in the nad/ directory of the proj source tree before proj is compiled and installed. From proj-datumgrids 1.6, the contents may also be added to share/proj, or maybe to sfr/proj, but historically grid *.lla files needed conversion with nad2bin during the installation of proj itself.

See current proj-datumgrid, with many other candidate grids that users can employ here.

For Windows, the files are copied in rgdal using rgdal/configure.win, and the same could be done for sfr. For OSX (as a Unix variant but built static), look at rgdal/configure.in lines from 415 and so on - they were written for Simon Urbanek. So configure args --with-data-copy=yes should trigger what configure.win does on Windows, that is adding the installed share/proj and share/gdal to R_PACKAGE_DIR.

Compiler warnings

On R 3.3 on Windows. The sign-compare warnings can be fixed by using size_t i instead of int i in for(int i = 0; ........ ) loops

"C:/PROGRA~1/R/R-33~1.1/bin/x64/Rscript.exe" "../tools/winlibs.R"
c:/Rtools/mingw_64/bin/g++  -std=c++0x -I"C:/PROGRA~1/R/R-33~1.1/include" -DNDEBUG -I../windows/gdal-2.0.1/include -I../windows/geos-3.5.0/include   -I"C:/Users/Jeroen/Documents/R/win-library/3.3/Rcpp/include" -I"d:/Compiler/gcc-4.9.3/local330/include"     -O2 -Wall  -mtune=core2 -c RcppExports.cpp -o RcppExports.o
c:/Rtools/mingw_64/bin/g++  -std=c++0x -I"C:/PROGRA~1/R/R-33~1.1/include" -DNDEBUG -I../windows/gdal-2.0.1/include -I../windows/geos-3.5.0/include   -I"C:/Users/Jeroen/Documents/R/win-library/3.3/Rcpp/include" -I"d:/Compiler/gcc-4.9.3/local330/include"     -O2 -Wall  -mtune=core2 -c bbox.cpp -o bbox.o
c:/Rtools/mingw_64/bin/g++  -std=c++0x -I"C:/PROGRA~1/R/R-33~1.1/include" -DNDEBUG -I../windows/gdal-2.0.1/include -I../windows/geos-3.5.0/include   -I"C:/Users/Jeroen/Documents/R/win-library/3.3/Rcpp/include" -I"d:/Compiler/gcc-4.9.3/local330/include"     -O2 -Wall  -mtune=core2 -c gdal.cpp -o gdal.o
gdal.cpp: In function 'Rcpp::List sfc_from_geometries(std::vector<OGRGeometry*>, bool)':
gdal.cpp:95:20: warning: comparison between signed and unsigned integer expressions [-Wsign-compare]
  for (int i = 0; i < g.size(); i++) {
                    ^
gdal.cpp: In function 'Rcpp::CharacterVector p4s_from_spatial_reference(OGRSpatialReference*)':
gdal.cpp:109:9: warning: unused variable 'err' [-Wunused-variable]
  OGRErr err = ref->exportToProj4(&cp);
         ^
gdal.cpp: In function 'Rcpp::List CPL_transform(Rcpp::List, Rcpp::CharacterVector)':
gdal.cpp:130:20: warning: comparison between signed and unsigned integer expressions [-Wsign-compare]
  for (int i = 0; i < g.size(); i++)
                    ^
c:/Rtools/mingw_64/bin/g++  -std=c++0x -I"C:/PROGRA~1/R/R-33~1.1/include" -DNDEBUG -I../windows/gdal-2.0.1/include -I../windows/geos-3.5.0/include   -I"C:/Users/Jeroen/Documents/R/win-library/3.3/Rcpp/include" -I"d:/Compiler/gcc-4.9.3/local330/include"     -O2 -Wall  -mtune=core2 -c gdal_read.cpp -o gdal_read.o
c:/Rtools/mingw_64/bin/g++  -std=c++0x -I"C:/PROGRA~1/R/R-33~1.1/include" -DNDEBUG -I../windows/gdal-2.0.1/include -I../windows/geos-3.5.0/include   -I"C:/Users/Jeroen/Documents/R/win-library/3.3/Rcpp/include" -I"d:/Compiler/gcc-4.9.3/local330/include"     -O2 -Wall  -mtune=core2 -c geos.cpp -o geos.o
In file included from ../windows/geos-3.5.0/include/geos/geom/Geometry.h:25:0,
                 from geos.cpp:6:
../windows/geos-3.5.0/include/geos/platform.h:96:0: warning: "ISNAN" redefined
 #  define ISNAN(x) (std::isnan(x))
 ^
In file included from C:/PROGRA~1/R/R-33~1.1/include/R.h:66:0,
                 from C:/Users/Jeroen/Documents/R/win-library/3.3/Rcpp/include/Rcpp/r/headers.h:52,
                 from C:/Users/Jeroen/Documents/R/win-library/3.3/Rcpp/include/RcppCommon.h:29,
                 from C:/Users/Jeroen/Documents/R/win-library/3.3/Rcpp/include/Rcpp.h:27,
                 from geos.cpp:4:
C:/PROGRA~1/R/R-33~1.1/include/R_ext/Arith.h:70:0: note: this is the location of the previous definition
 #  define ISNAN(x)     R_isnancpp(x)
 ^
In file included from ../windows/geos-3.5.0/include/geos/io/WKBReader.h:25:0,
                 from geos.cpp:7:
../windows/geos-3.5.0/include/geos/geom/GeometryFactory.h:209:49: warning: 'auto_ptr' is deprecated (declared at C:/Rtools/mingw_64/x86_64-w64-mingw32/include/c++/backward/auto_ptr.h:87) [-Wdeprecated-declarations]
    std::auto_ptr<CoordinateSequence> newCoords) const;
                                                 ^
../windows/geos-3.5.0/include/geos/geom/GeometryFactory.h:258:51: warning: 'auto_ptr' is deprecated (declared at C:/Rtools/mingw_64/x86_64-w64-mingw32/include/c++/backward/auto_ptr.h:87) [-Wdeprecated-declarations]
    std::auto_ptr<CoordinateSequence> coordinates) const;
                                                   ^
../windows/geos-3.5.0/include/geos/geom/GeometryFactory.h:306:58: warning: 'auto_ptr' is deprecated (declared at C:/Rtools/mingw_64/x86_64-w64-mingw32/include/c++/backward/auto_ptr.h:87) [-Wdeprecated-declarations]
  std::auto_ptr<Geometry> buildGeometry(T from, T toofar) const
                                                          ^
../windows/geos-3.5.0/include/geos/geom/GeometryFactory.h: In member function 'std::auto_ptr<geos::geom::Geometry> geos::geom::GeometryFactory::buildGeometry(T, T) const':
../windows/geos-3.5.0/include/geos/geom/GeometryFactory.h:325:66: warning: 'auto_ptr' is deprecated (declared at C:/Rtools/mingw_64/x86_64-w64-mingw32/include/c++/backward/auto_ptr.h:87) [-Wdeprecated-declarations]
       return std::auto_ptr<Geometry>( createGeometryCollection() );
                                                                  ^
../windows/geos-3.5.0/include/geos/geom/GeometryFactory.h:330:56: warning: 'auto_ptr' is deprecated (declared at C:/Rtools/mingw_64/x86_64-w64-mingw32/include/c++/backward/auto_ptr.h:87) [-Wdeprecated-declarations]
       return std::auto_ptr<Geometry>( (*from)->clone() );
                                                        ^
../windows/geos-3.5.0/include/geos/geom/GeometryFactory.h:348:75: warning: 'auto_ptr' is deprecated (declared at C:/Rtools/mingw_64/x86_64-w64-mingw32/include/c++/backward/auto_ptr.h:87) [-Wdeprecated-declarations]
       return std::auto_ptr<Geometry>( createGeometryCollection(fromGeoms) );
                                                                           ^
../windows/geos-3.5.0/include/geos/geom/GeometryFactory.h:353:69: warning: 'auto_ptr' is deprecated (declared at C:/Rtools/mingw_64/x86_64-w64-mingw32/include/c++/backward/auto_ptr.h:87) [-Wdeprecated-declarations]
       return std::auto_ptr<Geometry>( createMultiPolygon(fromGeoms) );
                                                                     ^
../windows/geos-3.5.0/include/geos/geom/GeometryFactory.h:355:72: warning: 'auto_ptr' is deprecated (declared at C:/Rtools/mingw_64/x86_64-w64-mingw32/include/c++/backward/auto_ptr.h:87) [-Wdeprecated-declarations]
       return std::auto_ptr<Geometry>( createMultiLineString(fromGeoms) );
                                                                        ^
../windows/geos-3.5.0/include/geos/geom/GeometryFactory.h:357:67: warning: 'auto_ptr' is deprecated (declared at C:/Rtools/mingw_64/x86_64-w64-mingw32/include/c++/backward/auto_ptr.h:87) [-Wdeprecated-declarations]
       return std::auto_ptr<Geometry>( createMultiPoint(fromGeoms) );
                                                                   ^
../windows/geos-3.5.0/include/geos/geom/GeometryFactory.h:361:36: warning: 'auto_ptr' is deprecated (declared at C:/Rtools/mingw_64/x86_64-w64-mingw32/include/c++/backward/auto_ptr.h:87) [-Wdeprecated-declarations]
     return std::auto_ptr<Geometry>(); // nullptr
                                    ^
geos.cpp: In function 'Rcpp::LogicalVector CPL_geos_is_valid(Rcpp::List)':
geos.cpp:139:11: warning: 'i' may be used uninitialized in this function [-Wmaybe-uninitialized]
  for (int i; i < sfc.length(); i++)
           ^
c:/Rtools/mingw_64/bin/g++  -std=c++0x -I"C:/PROGRA~1/R/R-33~1.1/include" -DNDEBUG -I../windows/gdal-2.0.1/include -I../windows/geos-3.5.0/include   -I"C:/Users/Jeroen/Documents/R/win-library/3.3/Rcpp/include" -I"d:/Compiler/gcc-4.9.3/local330/include"     -O2 -Wall  -mtune=core2 -c wkb.cpp -o wkb.o
wkb.cpp: In function 'Rcpp::NumericMatrix read_multipoint(unsigned char**, int, bool, int, Rcpp::CharacterVector)':
wkb.cpp:80:20: warning: comparison between signed and unsigned integer expressions [-Wsign-compare]
  for (int i = 0; i < *npts; i++) {
                    ^
wkb.cpp: In function 'Rcpp::List read_geometrycollection(unsigned char**, int, bool, int, Rcpp::CharacterVector, bool)':
wkb.cpp:96:20: warning: comparison between signed and unsigned integer expressions [-Wsign-compare]
  for (int i = 0; i < *nlst; i++)
                    ^
wkb.cpp: In function 'Rcpp::NumericMatrix read_numeric_matrix(unsigned char**, int, Rcpp::CharacterVector)':
wkb.cpp:120:17: warning: comparison between signed and unsigned integer expressions [-Wsign-compare]
  for (int i=0; i<(*npts); i++)
                 ^
wkb.cpp: In function 'Rcpp::List read_matrix_list(unsigned char**, int, Rcpp::CharacterVector)':
wkb.cpp:133:20: warning: comparison between signed and unsigned integer expressions [-Wsign-compare]
  for (int i = 0; i < (*nlst); i++)
                    ^
wkb.cpp: In function 'Rcpp::List read_data(unsigned char**, bool, int, bool, int*, uint32_t*)':
wkb.cpp:151:19: warning: unused variable 'dim' [-Wunused-variable]
  int sf_type = 0, dim = 0, n_dims = 0;
                   ^
wkb.cpp: In function 'void write_matrix(std::ostringstream&, Rcpp::NumericMatrix, double)':
wkb.cpp:372:29: warning: comparison between signed and unsigned integer expressions [-Wsign-compare]
  for (unsigned int i = 0; i < mat.nrow(); i++)
                             ^
wkb.cpp:373:30: warning: comparison between signed and unsigned integer expressions [-Wsign-compare]
   for (unsigned int j = 0; j < mat.ncol(); j++)
                              ^
wkb.cpp: In function 'Rcpp::List CPL_write_wkb(Rcpp::List, bool, int, Rcpp::CharacterVector, double)':
wkb.cpp:501:21: warning: comparison between signed and unsigned integer expressions [-Wsign-compare]
   for (int j = 0; j < str.size(); j++)
                     ^
wkb.cpp:487:6: warning: unused variable 'type' [-Wunused-variable]
  int type = 0, last_type = 0, n_types = 0;
      ^
wkb.cpp:487:16: warning: unused variable 'last_type' [-Wunused-variable]
  int type = 0, last_type = 0, n_types = 0;
                ^
wkb.cpp:487:31: warning: unused variable 'n_types' [-Wunused-variable]
  int type = 0, last_type = 0, n_types = 0;
                               ^

Encoding of layer names, attribute fields, text attributes

The links to OGC documents point to the geometries, but don't cover the question of the encoding representation of file names, layer names, attribute field names or text attributes. There are some hooks in OGR/CPL for this, which should be aligned with the cross-platform mechanisms in R. GeoJSON seems to want Unicode, GPKG maybe UTF-8, etc.

st_write_db enhancements

Maybe this function could use 'drop if exists' if available or wrap 'drop table' in a tryCatch.
Also, it currently doesn't seem to be able to append or update an existing table.

Here are a couple other things that I expected to work:

#does not pass schema/table on to dbWriteTable
st_write_db(conn = pg_con, obj = nc, table = c("public", "sf_test_nc"))

# does not pass row.names on to dbWriteTable
st_write_db(conn = pg_con, obj = nc, table = "sf_test_nc", row.names = FALSE)

Don't coerce to GEOMETRYCOLLECTION in case of mixed geometries

If we have more than one geometry type, do not coerce everything into a GEOMETRYCOLLECTION, but leave as is, and give the column type sfc_GEOMETRY. Allow a path to select particular geometry types. Add a promote_to_multi option to st_read (see ogr2ogr).

[Request] dplyr join operations with spatial predicates?

I haven't dug into sfr yet, so I don't know if you're allowing sf objects to contain multiple geometries (i.e. two MULTIPOLYGONs or a MULTIPOLYGON and a MULTIPOINT). If you can, (and if you can't, I think you should), then *join operations would result in sf objects with multiple geometries.

A proposed combo:

mp1 # sf with MULTIPOLYGON
mp2 # sf with MULTIPOLYGON
mp3 <- left_join(mp1, mp2, border = FALSE)

would result in an sf with two columns geometry1 and geometry2 and would behave just like left_join in dplyr except a match is defined as two geometries with positive spatial overlap; border = TRUE would include matches with common borders. left_join(x, y) includes all observations in x, regardless of whether they match or not. Multiple matches return as multiple rows.

There are maybe more match conditions that you would want to include. Perhaps allowing any function that takes two sfg objects and returns a logical would make the most sense. Does st_intersects, for example, allow sfg objects as inputs? If so, then this could be fully general for any pair of sf types. Usage would be something like:

mp3 <- left_join(mp1, mp2, match_fun = st_touches)

It probably also makes sense to do this so you don't have to reinvent the join wheel. I don't believe dplyr allows for arbitrary conditions for join matches yet (tidyverse/dplyr#557), but there might be a solution that isn't too difficult to implement.

The user could also generate new geometry columns using any function that takes two sfgs and returns an sfg. i.e.:

mp3 <- mp3 %>%
  mutate(geometry3 = purrr::map2(geometry1, geometry2, st_intersection))

Please test for minimum version libraries in configure

Congrats on getting sf onto CRAN.

With it being there, it got caught up in a reverse dependency check for Rcpp which promptly failed leading me to suspect we had broken something.

But we didn't. Your build system is too lenient. Your README.md makes it clear (via the pointer to the travis.yml) that certain minimal versions are required. Apparently my recent-ish Ubuntu machine with 16.04 is not good enough. But you can actually know what version your are using so I think you should use this and test for a minimal version. I can probably help you with some configure.in snippets if you need them.

Curve and surface in OGR SFS

Heads up for this ticket and associated thread

https://trac.osgeo.org/gdal/ticket/6401

This is for the Z/M/ZM changes, and related to getting GPKG up to full spec.

st_read example in the sf vignette

I've got an error when I tired to replicate results from the sf vignette:

library('sf')
u_kmz <- "http://coagisweb.cabq.gov/datadownload/BikePaths.kmz"
download.file(u_kmz, "BikePaths.kmz")
biketrails_kmz <- st_read("BikePaths.kmz")

Cannot open data source BikePaths.kmz
Error in eval(substitute(expr), envir, enclos) : Open failed.

(based on today's version of the sf package from github)

scales?

I have heard through the grapevine (I'm still on a twitter sabbatical) that there's need of a ggplot scale for something geo-related (I'm assuming geo-related). Is this something @thomasp85 or I can help with?

empty geometries

GDAL seems to return empty geoms as GEOMETRYCOLLECTION with zero elements:

> st_difference(st_sfc(st_point(0:1),st_point(c(0,0))), st_sfc(st_point(0:1)))
Geometry set for 2 features 
geometry type:  GEOMETRYCOLLECTION
dimension:      XY
bbox:           xmin: 0 ymin: 0 xmax: 0 ymax: 0
epsg (SRID):    NA
proj4string:    NA
precision:      double (default; no precision model)
GEOMETRYCOLLECTION() 
GEOMETRYCOLLECTION(POINT(0 0))

The current logic of sfr is then to convert everything back to GEOMETRYCOLLECTION, only because the empty one is. Having empty points doesn't seem a good solution, as they don't have a WKB representation. Maybe tweak sfc such that they have a constant type (in this case POINT), but may contain empty geometries of type GEOMETRYCOLLECTION?

Extract centroids coordinates as a matrix

What would be the correct way to extract centroids coordinates of sf objects (I have sp::coordinates in mind) ? There must something easier/cleaner than:

nc <- st_read(system.file("shape/nc.shp", package="sf"), "nc", crs = 4267)
st_geometry(nc) <- st_centroid(nc)
pts <- st_as_text(st_geometry(nc))
m <- matrix(data = as.numeric(unlist(strsplit(substr(pts, 7, nchar(pts)-1), split = " "))), ncol = 2, byrow = T)

path.expand dsn arguments?

Currently the dsn argument for st_read, st_write, and st_list are directly passed into gdal for processing. The issue is that gdal does not understand common bash directory shortcuts like ~/shapefile/ which results in confusing errors like:

Cannot open data source ~/shapefile/
Error in eval(substitute(expr), envir, enclos) : Open failed.

This is problematic for new users because they can see that clearly ~/shapefile/ does exist.

Is there any harm in making sure path.expand is used on each dsn argument before passing it along to gdal?

rgeos:readWKT study case

I am reporting here on Roger Bivand's advise about a strange behaviour with rgeos:readWKT.

Please consider the attached tmp1 variable (just extracted from QGIS attribute table by copy and paste and read.delim from the clipboard). First unzip tmp1.zip

load("tmp1.Rdata")
library(rgeos)
readWKT(tmp1)
readWKT(as.character(tmp1))
readWKT(gsub(" ","",as.character(tmp1)))

In every case I get:

Error: Unable to parse: Polygon ((987298.45400000002700835 6684766.64900000020861626, 987354.71499999996740371 6684746.2369999997317791, 987430.25699999998323619 6684717.31500000040978193, 987477.16800000006332994 6684712.3030000003054738, 987506.0400000000372529 6684709.22500000055879354, 987529.53000000002793968 6684706.72400000039488077, 987571.87600000004749745 6684731.8880000002682209, 987529.01000000000931323 6684831.39400000032037497, 987503.28200000000651926 6684914.84300000034272671, 987479.98499999998603016 6684960.95500000007450581, 987479.231000000028871 6684966.91200000047683716, 987473.37300000002142042 6684977.03500000014901161, 987466.03700000001117587 6684989.57899999991059303, 987531.43599999998696148 6684992.7900000000372529, 987563.93099999998230487 6685019.06400000024586916, 987634.05099999997764826 6685050.60900000017136335, 987679.39800000004470348 6684986.36400000005960464, 987723.43099999998230487 6685007.84599999990314245, 987765.4599999999627471 6685035.6409

Back to the doc, I found that a possible explanation is that as.character breaks lines in language objects at 500 characters, and inserts newlines.

In our case, I suppose factors may be coerced to character using as.character, thus will not permit WKT description > 500 characters. Now the question is to get a workaround: many spatial objects are described with more than 500 characters in WKT...

tmp1.zip

installation error: linked libraries not found in Windows

This may be my environment that's not set up correctly: it looks like a some libraries are not found when trying to install the package in Windows (64-bit). Here's the console output from RStudio:

Downloading GitHub repo edzer/sfr@master
from URL https://api.github.com/repos/edzer/sfr/zipball/master
Installing sf
Skipping 2 packages ahead of CRAN: DBI, Rcpp
"C:/PROGRA~1/MIE74D~1/MRO/R-32~1.4/bin/x64/R" --no-site-file --no-environ --no-save  \
  --no-restore CMD INSTALL  \
  "C:/Users/frogu/AppData/Local/Temp/RtmpSW6DYg/devtoolsfec4c1d2a14/edzer-sfr-e2495a0"  \
  --library="E:/Data/R/library" --install-tests 

* installing *source* package 'sf' ...
Loading required package: methods
Loading required package: Rcpp
Loading required package: DBI
in method for 'coerce' with signature '"sf","Spatial"': no definition for class "Spatial"
in method for 'coerce' with signature '"sfc","Spatial"': no definition for class "Spatial"
Writing st_bbox.Rd
Writing st_read.Rd
Writing st_write.Rd
Writing st_write_db.Rd
Writing sfc.Rd
Writing summary.sfc.Rd
Writing st_geometry.Rd
Writing tibble.Rd
Writing st.Rd
Writing st_as_sf.Rd
Writing sf.Rd
Writing st_as_sfc.Rd
Writing st_as_binary.Rd
Writing st_as_text.Rd
Writing plot.Rd
Writing geos.Rd
Writing st_transform.Rd
Writing crs.Rd
** libs
Warning: this package has a non-empty 'configure.win' file,
so building only the main architecture

rm -f RcppExports.o bbox.o gdal.o gdal_geom.o gdal_read.o gdal_write.o geos.o wkb.o sf.dll
mkdir -p ../inst
"C:/PROGRA~1/MIE74D~1/MRO/R-32~1.4/bin/x64/Rscript.exe" "../tools/winlibs.R"
cp -r ../windows/gdal2-2.1.1/share/gdal ../inst/
cp -r ../windows/gdal2-2.1.1/share/proj ../inst/
g++ -m64 -std=c++0x -I"C:/PROGRA~1/MIE74D~1/MRO/R-32~1.4/include" -DNDEBUG -I../windows/gdal2-2.1.1/include/gdal -I../windows/gdal2-2.1.1/include/geos   -I"E:/Data/R/library/Rcpp/include" -I"c:/applications/extsoft/include"     -O2 -Wall  -mtune=core2 -c RcppExports.cpp -o RcppExports.o
g++ -m64 -std=c++0x -I"C:/PROGRA~1/MIE74D~1/MRO/R-32~1.4/include" -DNDEBUG -I../windows/gdal2-2.1.1/include/gdal -I../windows/gdal2-2.1.1/include/geos   -I"E:/Data/R/library/Rcpp/include" -I"c:/applications/extsoft/include"     -O2 -Wall  -mtune=core2 -c bbox.cpp -o bbox.o
g++ -m64 -std=c++0x -I"C:/PROGRA~1/MIE74D~1/MRO/R-32~1.4/include" -DNDEBUG -I../windows/gdal2-2.1.1/include/gdal -I../windows/gdal2-2.1.1/include/geos   -I"E:/Data/R/library/Rcpp/include" -I"c:/applications/extsoft/include"     -O2 -Wall  -mtune=core2 -c gdal.cpp -o gdal.o
g++ -m64 -std=c++0x -I"C:/PROGRA~1/MIE74D~1/MRO/R-32~1.4/include" -DNDEBUG -I../windows/gdal2-2.1.1/include/gdal -I../windows/gdal2-2.1.1/include/geos   -I"E:/Data/R/library/Rcpp/include" -I"c:/applications/extsoft/include"     -O2 -Wall  -mtune=core2 -c gdal_geom.cpp -o gdal_geom.o
g++ -m64 -std=c++0x -I"C:/PROGRA~1/MIE74D~1/MRO/R-32~1.4/include" -DNDEBUG -I../windows/gdal2-2.1.1/include/gdal -I../windows/gdal2-2.1.1/include/geos   -I"E:/Data/R/library/Rcpp/include" -I"c:/applications/extsoft/include"     -O2 -Wall  -mtune=core2 -c gdal_read.cpp -o gdal_read.o
gdal_read.cpp: In function 'Rcpp::List CPL_read_ogr(Rcpp::CharacterVector, Rcpp::CharacterVector, Rcpp::CharacterVector, bool, int, int)':
gdal_read.cpp:186:35: warning: comparison between signed and unsigned integer expressions [-Wsign-compare]
g++ -m64 -std=c++0x -I"C:/PROGRA~1/MIE74D~1/MRO/R-32~1.4/include" -DNDEBUG -I../windows/gdal2-2.1.1/include/gdal -I../windows/gdal2-2.1.1/include/geos   -I"E:/Data/R/library/Rcpp/include" -I"c:/applications/extsoft/include"     -O2 -Wall  -mtune=core2 -c gdal_write.cpp -o gdal_write.o
gdal_write.cpp: In function 'void SetFields(OGRFeature*, std::vector<OGRFieldType>, Rcpp::List, size_t)':
gdal_write.cpp:51:37: warning: comparison between signed and unsigned integer expressions [-Wsign-compare]
g++ -m64 -std=c++0x -I"C:/PROGRA~1/MIE74D~1/MRO/R-32~1.4/include" -DNDEBUG -I../windows/gdal2-2.1.1/include/gdal -I../windows/gdal2-2.1.1/include/geos   -I"E:/Data/R/library/Rcpp/include" -I"c:/applications/extsoft/include"     -O2 -Wall  -mtune=core2 -c geos.cpp -o geos.o
In file included from ../windows/gdal2-2.1.1/include/geos/geos/geom/Geometry.h:25:0,
                 from geos.cpp:6:
../windows/gdal2-2.1.1/include/geos/geos/platform.h:96:0: warning: "ISNAN" redefined [enabled by default]
C:/PROGRA~1/MIE74D~1/MRO/R-32~1.4/include/R_ext/Arith.h:81:0: note: this is the location of the previous definition
In file included from ../windows/gdal2-2.1.1/include/geos/geos/io/WKBReader.h:25:0,
                 from geos.cpp:7:
../windows/gdal2-2.1.1/include/geos/geos/geom/GeometryFactory.h:209:49: warning: 'auto_ptr' is deprecated (declared at c:\mrnmicro\applic\rtools\gcc-4.6.3\bin\../lib/gcc/i686-w64-mingw32/4.6.3/../../../../include/c++/4.6.3/backward/auto_ptr.h:87) [-Wdeprecated-declarations]
../windows/gdal2-2.1.1/include/geos/geos/geom/GeometryFactory.h:258:51: warning: 'auto_ptr' is deprecated (declared at c:\mrnmicro\applic\rtools\gcc-4.6.3\bin\../lib/gcc/i686-w64-mingw32/4.6.3/../../../../include/c++/4.6.3/backward/auto_ptr.h:87) [-Wdeprecated-declarations]
../windows/gdal2-2.1.1/include/geos/geos/geom/GeometryFactory.h:306:58: warning: 'auto_ptr' is deprecated (declared at c:\mrnmicro\applic\rtools\gcc-4.6.3\bin\../lib/gcc/i686-w64-mingw32/4.6.3/../../../../include/c++/4.6.3/backward/auto_ptr.h:87) [-Wdeprecated-declarations]
g++ -m64 -std=c++0x -I"C:/PROGRA~1/MIE74D~1/MRO/R-32~1.4/include" -DNDEBUG -I../windows/gdal2-2.1.1/include/gdal -I../windows/gdal2-2.1.1/include/geos   -I"E:/Data/R/library/Rcpp/include" -I"c:/applications/extsoft/include"     -O2 -Wall  -mtune=core2 -c wkb.cpp -o wkb.o
g++ -m64 -shared -s -static-libgcc -o sf.dll tmp.def RcppExports.o bbox.o gdal.o gdal_geom.o gdal_read.o gdal_write.o geos.o wkb.o -L../windows/gdal2-2.1.1/lib-4.6.3/x64 -lgdal -lsqlite3 -lspatialite -lproj -lgeos_c -lgeos -lhdf5 -lexpat -lfreexl -lcfitsio -lpng16 -ljpeg -ltiff -lgeotiff -lgif -lxml2 -llzma -lszip -lodbc32 -lodbccp32 -liconv -lz -lws2_32 -Lc:/applications/extsoft/lib/x64 -Lc:/applications/extsoft/lib -LC:/PROGRA~1/MIE74D~1/MRO/R-32~1.4/bin/x64 -lR
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -lgdal
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -lsqlite3
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -lspatialite
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -lproj
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -lgeos_c
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -lgeos
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -lhdf5
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -lexpat
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -lfreexl
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -lcfitsio
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -lpng16
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -ljpeg
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -ltiff
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -lgeotiff
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -lgif
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -lxml2
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -llzma
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -lszip
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -liconv
c:/mrnmicro/applic/rtools/gcc-4.6.3/bin/../lib/gcc/i686-w64-mingw32/4.6.3/../../../../i686-w64-mingw32/bin/ld.exe: cannot find -lz
collect2: ld returned 1 exit status
aucune DLL n'a pas été créé
ERROR: compilation failed for package 'sf'
* removing 'E:/Data/R/library/sf'
Error: Command failed (1)

I tried installing both from RStudio and from the command line, and also as system admin. Here are more details of my setup:

> R.version
platform       x86_64-w64-mingw32          
arch           x86_64                      
os             mingw32                     
system         x86_64, mingw32             
status                                     
major          3                           
minor          2.4                         
year           2016                        
month          03                          
day            16                          
svn rev        70336                       
language       R                           
version.string R version 3.2.4 (2016-03-16)
nickname       Very Secure Dishes

Do you have an idea what might be missing or wrongly set in my setup/environment?

Thanks for helping me sorting this out.

[Request] dplyr::bind_rows() for sf

Building on the list of dplyr verbs requested in edzer/sfr#42, could you please consider adding dplyr::bind_rows()?

This function makes it easy to combine dataframe-like objects, which would also be useful for sf objects.

Unlike rbind() or spRbind(), the bind_rows() function allows the merger of objects with non-matching columns, filling any unshared columns with NA. I find this convenience feature saves me a lot of time, even if it does lead to the creation of the occasional ugly dataframe.

There are probably some details that would need to be worked out with the geom list-columns, especially in cases where sf objects with differing geometry classes are present. Perhaps you could coerce the geom column to geometry type: GEOMETRY?

empty raw WKB policy

With an empty WKB element I get errors, whether by pureR or not:

library(RODBC)
## not reproducible without Manifold and mdsumner/manifoldr installed
mfile <- "... ManifoldDVD\\Data\\Europe\\Hydrography.map"
con <- manifoldr::odbcConnectManifold(mfile)
x <- sqlQuery(con, "SELECT ID, CGeomWKB(Geom(ID)) AS [geom]  FROM [Hydrography Table]")
library(sf)
sf_geom_name <- "geom"
 
x[[sf_geom_name]] <- sf::st_as_sfc(structure(x[[sf_geom_name]], class = "WKB"), EWKB = FALSE)

#Error in eval(substitute(expr), envir, enclos) : 
 # non native endian: use pureR = TRUE

x[[sf_geom_name]] <- sf::st_as_sfc(structure(x[[sf_geom_name]], class = "WKB"), EWKB = FALSE, pureR = T)
#Error in readBin(rc, what = "integer", n = 1L, size = 4L, endian = endian) : 
#  invalid 'swap' argument

Here's what the empty geom and its neighbour look like:

x$geom[55285:55286]
[[1]]
raw(0)

[[2]]
  [1] 01 03 00 00 00 01 00 00 00 09 00 00 00 00 00 00 00 bc 43 10 41 00 00 00 00 c6 33 23 41 00 00 00 00 cc 43 10 41 00 00 00 00 7a 33 23 41 00 00
 [48] 00 00 38 46 10 41 00 00 00 00 7a 32 23 41 00 00 00 00 6c 49 10 41 00 00 00 00 a2 31 23 41 00 00 00 00 60 4a 10 41 00 00 00 00 c4 31 23 41 00
 [95] 00 00 00 78 4a 10 41 00 00 00 00 12 32 23 41 00 00 00 00 bc 48 10 41 00 00 00 00 ce 32 23 41 00 00 00 00 68 46 10 41 00 00 00 00 16 33 23 41
[142] 00 00 00 00 bc 43 10 41 00 00 00 00 c6 33 23 41

Should sf handle these broken cases? Is there any policy like this in the postgres work already for NULL geometries?

I'll pursue the Manifold cast to WKB as a external problem, but it seems this could happen more generally from DB sources.

I tried to make a reproducible example, independent of Manifold but can't get it to work, sorry.

Feature request -- eliding and ± scaling multipolygons; PLUS need more than just *.shp

First, the minor gotcha, for the benefit of other pilgrims. I was trying to read in a shapefile

filename <- system.file("/Users/rc/Desktop/gitmo/guardian/data/cb_2015_us_cd114_5m.shp", package="sf")

and got

   Error in default_layer(dsn) : 
  cannot guess layer name from this datasource

but then realized my data was not readable with system.file because it wasn't in a package.

Then I tried

    filename <- "/Users/rc/Desktop/gitmo/guardian/data/cb_2015_us_cd114_5m.shp"
    > us <- st_read(filename, crs = 4267)
      Cannot open data source /Users/rc/Desktop/gitmo/guardian/data/cb_2015_us_cd114_5m.shp
      Error in eval(substitute(expr), envir, enclos) : Open failed.
    In addition: Warning message:
      In CPL_read_ogr(dsn, layer, as.character(options), quiet, iGeomField -  :
     GDAL Error 4: Unable to open /Users/rc/Desktop/gitmo/guardian/data/cb_2015_us_cd114_5m.shx             or /Users/rc/Desktop/gitmo/guardian/data/cb_2015_us_cd114_5m.SHX.Try --config      SHAPE_RESTORE_SHX true

Finally, I realized that .shp still needs to be in the same directory with the other ESRI files

> filename <- "/Users/rc/Desktop/gitmo/guardian/data/cb_2015_us_cd114_5m/cb_2015_us_cd114_5m.shp"
> us <- st_read(filename, crs = 4267)
Reading layer cb_2015_us_cd114_5m from data source /Users/rc/Desktop/gitmo/guardian/data/cb_2015_us_cd114_5m/cb_2015_us_cd114_5m.shp using driver "ESRI Shapefile"
features:       441
fields:         8
converted into: MULTIPOLYGON
proj4string:    +proj=longlat +datum=NAD83 +no_defs 
Warning message:
In check_replace(x, value) :
  st_crs: replacing crs does not reproject data; use st_transform for that

So, in this way, it works just like ReadOGR.

The feature I need are the ability to elide and scale (as shown in the US map in the June presentation (not obvious how to do this in sf) and to be able to scale a multipolygon up. My objective is to outset the really small units like RI, DE and especially DC out in the Atlantic, big enough to add centroid labels and to color with attributes.

r-spatial / sf Goto Github PK

sf's Introduction

Simple Features for R

Books, journal articles, blogs, presentations, vignettes, sp-sf wiki

Cheatsheet

Installing

Windows

MacOS

Linux

Ubuntu

Fedora

Arch

renv or conda

Other

Multiple GDAL, GEOS and/or PROJ versions on your system

lwgeom

Contributing

How to cite

Acknowledgment

sf's People

Contributors

Stargazers

Watchers

Forkers

sf's Issues

Recommend Projects

Recommend Topics

Recommend Org

`renv` or `conda`