3

我想执行一个非常简单的操作：合并两个形状文件。具体来说，我有美国每个州的人口普查区形状文件，我想将它们合并到一个形状文件中。最终，我想获取组合的形状文件并在一组纬度和经度坐标上执行叠加，以确定我的坐标落入哪些人口普查区域。

我已经看到了很多关于这个的讨论（在 R 中结合边界 shapefile）。然而，所有的讨论都已经过时了，我希望在此期间这些软件包得到了改进。

我正在使用的文件位于此处： ftp: //ftp2.census.gov/geo/tiger/TIGER2010/TRACT/2010/

下面的代码可用于重新创建我正在使用的文件。它需要下载两个文件，总共大约 11 兆字节。运行时间应该只有一分钟。

任何帮助是极大的赞赏。这似乎是一个微不足道的操作。也许如果我在地理空间映射方面有更多经验，我可以更好地利用可用的文档。

以下是我尝试过的几件事：

### Insert your file path here
FPATH <- './data'

### Set up library
require(rgeos)
require(maptools)
require(RCurl)
require(parallel)
cl <- makeCluster(detectCores())

### Download files... (~11,000 KB total for this example)
ftp <- 'ftp://ftp2.census.gov/geo/tiger/TIGER2010/TRACT/2010/'
files <- getURLContent(ftp, dirlistonly = T)
files <- unlist(strsplit(files, split = '\r\n', fixed = T))
files <- grep('2010_[[:digit:]]{2}_', files, value = T)[1:2]  # Only grab two files for this example
clusterMap(cl, download.file, url = paste0(ftp, files), destfile = paste0(FPATH, files))

### Unzip shape files
files <- list.files(FPATH, full.names = T)
clusterMap(cl, unzip, zipfile = files, exdir = FPATH)

### Read in shape files
files <- list.files(FPATH, pattern = "shp$", full.names = T)
a <- readShapePoly(fn = files[1])
b <- readShapePoly(fn = files[2])

### Attempt to join two shape files
spRbind(a, b)  # Error in spRbind(as(obj, "SpatialPolygons"), as(x, "SpatialPolygons")) : non-unique polygon IDs
gUnion(a, b)   # Error in RGEOSBinTopoFunc(spgeom1, spgeom2, byid, id, "rgeos_union")   : geos_geospolygon2SpatialPolygons: ng > length(IDs)

感谢您的时间。

4

1 回答 1

5

You've got duplicated polygon IDs. You can change them manually, or use taRifx.geo:

library(devtools)
install_git("git://github.com/gsk3/taRifx.geo.git")
library(taRifx.geo)
rbind(a,b, fix.duplicated.IDs=TRUE)

Code available for inspection here.

于 2013-11-13T19:26:39.650 回答