时间戳

首页 > 解决方案 > 使用 sf 包在 R 中读取大型 shapefile

r - 使用 sf 包在 R 中读取大型 shapefile

问题描述

我有一个带有 32M 观察点的 shapefile。我想将它加载到 R 中，我尝试了 read_sf 和 st_read 但我的 R 会话不断崩溃。我想到的另一种方法是编写一个 for 循环，对我想要的列进行子集化，并且可能一次执行特定数量的行，然后对它们进行 rbinding，但无法弄清楚如何让 R 理解查询。这是我到目前为止不起作用的内容：

for (i in 1:10) {
  j = i-1
  jj = i+1
  print(i)
  print(j)
  print(jj)
  A <- read_sf("C:\\Users\\...parcels-20210802T125336Z-001\\parcels\\join_L3_Mad_Addresses.shp", query = "SELECT FID, CENTROID_I, LOC_ID FROM join_L3_Mad_Addresses WHERE FID < "jj" AND FID > "j"")
}

标签： rshapefilesf

解决方案

我认为您可以重新编写以下代码。

加载包

library(sf)
#> Linking to GEOS 3.9.0, GDAL 3.2.1, PROJ 7.2.1

定义 shapefile 的路径

dsn <- system.file("shape/nc.shp", package="sf")

计算特征的数量dsn

st_layers(dsn, do_count = TRUE)
#> Driver: ESRI Shapefile 
#> Available layers:
#>   layer_name geometry_type features fields
#> 1         nc       Polygon      100     14

启动一个 for 循环，一次读取 10 个特征。将数据添加到列表

shp_data_list <- list()
i <- 1
for (offset in seq(10, 100, by = 10)) {
  query <- paste0("SELECT * FROM nc LIMIT ", 10, " OFFSET ", offset - 10)
  shp_data_list[[i]] <- st_read(dsn, query = query, quiet = TRUE)
  gc(verbose = FALSE)
  i <- i + 1
}

绑定对象

shp_data <- do.call(rbind, shp_data_list)

添加一个 ID 列（仅用于绘图）

shp_data$ID <- as.character(rep(1:10, each = 10))
plot(shp_data["ID"])

唯一的问题是这个过程可能不会保留几何类型。例如，

unique(st_geometry_type(shp_data))
#> [1] MULTIPOLYGON POLYGON     
#> 18 Levels: GEOMETRY POINT LINESTRING POLYGON MULTIPOINT ... TRIANGLE

尽管

unique(st_geometry_type(st_read(dsn, quiet = TRUE)))
#> [1] MULTIPOLYGON
#> 18 Levels: GEOMETRY POINT LINESTRING POLYGON MULTIPOINT ... TRIANGLE

您可以更改几何类型st_cast()

由reprex 包于 2021-08-03 创建 (v2.0.0 )

推荐阅读