将“多特征”GeoJSON 的特征转换为 R 空间对象答案

【问题标题】：Convert features of a 'multifeature' GeoJSON into R spatial objects将“多特征”GeoJSON 的特征转换为 R 空间对象
【发布时间】：2015-06-01 20:50:42
【问题描述】：

通常您可以使用可靠的readOGR 将 geojson 文件读入 R，如图所示 here。

但是，对于多特征 geojson，这会失败。

可重现的例子：

downloader::download("https://github.com/Robinlovelace/Creating-maps-in-R/raw/master/data/test-multifeature.geojson", "test.geojson")
test <- rgdal::readOGR("test.geojson", "OGRGeoJSON") # fails with:

Error in ogrInfo(dsn = dsn, layer = layer, encoding = encoding, use_iconv = use_iconv,  : 
  Multiple incompatible geometries: wkbPoint: 98; wkbLineString: 660; wkbPolygon: 23

错误信息很清楚，并指出了解决方案：拆分特征。但是，除了使用正则表达式之外，我不知道该怎么做。

非常欢迎任何想法。

令人惊奇的是：GitHub displays the data natively on the browser，而 R 甚至无法（似乎）读入它！

解决方案的替代方法：

test <- geojsonio::geojson_read("test.geojson")

【问题讨论】：

标签： r geojson gdal rgdal ropensci

【解决方案1】：

您可以对各种GDAL函数使用require_geomType参数来提取您需要的特征：

library(rgdal)

ogrListLayers("test.geojson")
## [1] "OGRGeoJSON"
## attr(,"driver")
## [1] "GeoJSON"
## attr(,"nlayers")
## [1] 1

# This fails but you can at least see the geoms it whines about
ogrInfo("test.geojson", "OGRGeoJSON")
## Error in ogrInfo("test.geojson", "OGRGeoJSON") : 
##   Multiple incompatible geometries: wkbPoint: 98; wkbLineString: 660; wkbPolygon: 23


ogrInfo("test.geojson", "OGRGeoJSON", require_geomType="wkbPoint")
## NOTE: keeping only 98 wkbPoint of 781 features
##     note that extent applies to all features
## Source: "test.geojson", layer: "OGRGeoJSON"
## Driver: GeoJSON number of rows 781 
##   selected geometry type: wkbPoint with 98 rows
## Feature type: wkbPoint with 2 dimensions
## Extent: (12.48326 41.88355) - (12.5033 41.89629)
## CRS: +proj=longlat +datum=WGS84 +no_defs  
## Number of fields: 78 
##                        name type length typeName
## 1                      area    4      0   String
## 2                   bicycle    4      0   String
## ...
## LONG LIST - 78 total


ogrInfo("test.geojson", "OGRGeoJSON", require_geomType="wkbLineString")
## NOTE: keeping only 660 wkbLineString of 781 features
##     note that extent applies to all features
## Source: "test.geojson", layer: "OGRGeoJSON"
## Driver: GeoJSON number of rows 781 
##   selected geometry type: wkbLineString with 660 rows
## Feature type: wkbLineString with 2 dimensions
## Extent: (12.48326 41.88355) - (12.5033 41.89629)
## CRS: +proj=longlat +datum=WGS84 +no_defs  
## Number of fields: 78 
##                        name type length typeName
## 1                      area    4      0   String
## 2                   bicycle    4      0   String
## ...
## LONG LIST - 78 total (same as above)


ogrInfo("test.geojson", "OGRGeoJSON", require_geomType="wkbPolygon")
## NOTE: keeping only 23 wkbPolygon of 781 features
##     note that extent applies to all features
## Source: "test.geojson", layer: "OGRGeoJSON"
## Driver: GeoJSON number of rows 781 
##   selected geometry type: wkbPolygon with 23 rows
## Feature type: wkbPolygon with 2 dimensions
## Extent: (12.48326 41.88355) - (12.5033 41.89629)
## CRS: +proj=longlat +datum=WGS84 +no_defs  
## Number of fields: 78 
##                        name type length typeName
## 1                      area    4      0   String
## 2                   bicycle    4      0   String
## ...
## LONG LIST - 78 total (same as above)


points <- readOGR("test.geojson", "OGRGeoJSON", require_geomType="wkbPoint")
## OGR data source with driver: GeoJSON 
## Source: "test.geojson", layer: "OGRGeoJSON"
## with 781 features;
## Selected wkbPoint feature type, with 98 rows
## It has 78 fields
## NOTE: keeping only 98 wkbPoint of 781 features

lines <- readOGR("test.geojson", "OGRGeoJSON", require_geomType="wkbLineString")
## OGR data source with driver: GeoJSON 
## Source: "test.geojson", layer: "OGRGeoJSON"
## with 781 features;
## Selected wkbLineString feature type, with 660 rows
## It has 78 fields
## NOTE: keeping only 660 wkbLineString of 781 features

polygons <- readOGR("test.geojson", "OGRGeoJSON", require_geomType="wkbPolygon")
## OGR data source with driver: GeoJSON 
## Source: "test.geojson", layer: "OGRGeoJSON"
## with 781 features;
## Selected wkbPolygon feature type, with 23 rows
## It has 78 fields
## NOTE: keeping only 23 wkbPolygon of 781 features

# prove they red in things
plot(lines, col="#7f7f7f")
plot(polygons, add=TRUE)
plot(points, add=TRUE, col="red")

【讨论】：

这是 rgdal 的最新版本吗？我在 0.9-1 中没有看到 - 刚刚升级到 0.9-3，我现在拥有它！
我认为它在 0.9-2 中。它非常有用，尤其是对于这些天酷孩子们正在使用的所有这些疯狂的 GeoJSON shapefile ;-)
酷。可以在这里找到我对 R/gdal 的多特征/GeoJSON 探索的总结：rpubs.com/RobinLovelace/84577。如果你使用 Ubuntu，这可以更新 gdal： $ sudo add-apt-repository ppa:ubuntugis/ubuntugis-unstable && sudo apt-get update $ sudo apt-get install gdal-bin

【解决方案2】：

你可以在命令行中使用 ogr2ogr 将这个巨大的嵌合体拆分成合理的东西：

ogr2ogr -where "OGR_GEOMETRY='LINESTRING'" \
     -f "GeoJSON" lines.geojson  test.geojson

同样适用于点和多边形。

几年前有一些关于在readOGR 中实现 OGR_SQL 的讨论，此时您可以从 R 中执行此操作，但 Roger 不想这样做，也没有人愿意提供帮助 :(

创建拆分的 geojson 文件后，您可以将它们读入单个 rgeos::SpatialCollections 对象：

points=readOGR("points.geojson","OGRGeoJSON")
polys=readOGR("polygons.geojson","OGRGeoJSON")
lines=readOGR("lines.geojson","OGRGeoJSON")
require(rgeos)
g = SpatialCollections(points=points, lines=lines, polygons=polys)
plot(g)

如果您想尝试使用geojsonio，那么您可以使用Filter 从几何集合中选择给定几何的列表元素

polygon_features = Filter(
    function(f){f$geometry$type=="Polygon"},
    test$features)

但是你仍然必须构建一些可以进入单独的 R 实体的东西....

【讨论】：

【解决方案3】：

几年后，一些替代方案 - library(geojsonsf) 和 library(sf) 都将读取 geojson 并转换为 sf 对象

url <- 'https://github.com/Robinlovelace/Creating-maps-in-R/raw/master/data/test-multifeature.geojson'

## these give the same result
sf <- geojsonsf::geojson_sf( url )
sf <- sf::st_read( url )

我们来看看

library(mapdeck)

set_token( "MAPBOX_TOKEN" )

mapdeck( style = mapdeck_style("light") ) %>%
    add_sf( sf )

【讨论】：