與ESRI形狀文件

多邊形匹配GIS點I具有含約150個連續的地理區域（多邊形），它們一起構成一個地理區域的庫ESRI形狀文件。我也有一個csv文件，其中包含60,000個事件，每個事件都有一組唯一的x，y點座標。在形狀文件中的150個多邊形的一個（也是唯一一個）中的每一個CSV文件的事件發生，但我不知道該多邊形的一個與每個記錄相關聯。因此，我需要編寫一個算法，以找出其中每個事件發生的多邊形的身份。我需要編寫的算法輸出將使我能夠隨後計算統計數據，例如特定多邊形（地理區域）內發生特定類型事件的可能性。（當然，形狀文件不僅僅是一個文件，它是一個包含8個文件的目錄，其文件擴展名包括.dbf，.prj，.qix，.sbn，.sbx，.shp，.xml和.SHX。）與ESRI形狀文件

我沒有一個ArcGIS許可。我已經找到了文件地理數據庫API在http://resources.arcgis.com/content/geodatabases/10.0/file-gdb-api，但我不知道它是正確的工具集，而我也無法找到示例代碼。

任何人都可以看到一些用於查找哪個多邊形（來自形狀文件）的大量點（來自諸如csv文件之類的外部數據源）的代碼落入內部嗎？

此外，我需要指定我需要能夠爲每個事件的csv記錄添加相關多邊形的標識的特定代碼。因此，僅僅繪製地圖上的點來顯示哪些多邊形包含事件是不夠的。我根本不需要可視化這些數據。相反，我需要的是能夠標記多邊形ID每個事件記錄在CSV文件中，這樣我可以做後續的數值分析，是不是在視覺性。

文章鏈接，教程，並在此主題的其他資源也很受青睞。我認爲這是一個人們每天都在解決的問題，因此如果一個人知道如何去尋找它，那麼必須有一個確定的代碼庫。我每天都用Java編寫代碼，因此Java解決方案是首選。但是，如果您有使用不同語言編寫的良好代碼示例，則可以從另一種語言移植某些內容。

* 編輯：*
我試圖根據Spacedman的建議下列R-代碼，我得到了以下錯誤消息：

> myCSV <- read.csv(file="myCSVFile.csv",head=TRUE,sep=",") 
> pts = SpatialPoints(myCSV) 
> ZipCodes = readShapeSpatial("path/myshapefile.shp") 
> overlay(myCSV,ZipCodes) 
Error in function (classes, fdef, mtable) : unable to find an inherited method for function "overlay", for signature "data.frame", "SpatialPolygonsDataFrame" 
>

見下面我的其他意見。

第二個編輯：
，我結束了使用R代碼爲：

myCSV <- read.csv(file="myData.csv",head=TRUE,sep=",") 
pts = SpatialPoints(myCSV) 
ZipCodes = readShapeSpatial("myPath/ZipCodes.shp") 
write.csv(ZipCodes$ZIPCODE[overlay(pts,ZipCodes)], "ZipMatches.csv", quote=FALSE, row.names=FALSE)

注：我不得不使用：

summary(ZipCodes)

找到包含郵政編碼編碼的字段的名稱。在我運行摘要（ZipCodes）之前，腳本只是放出每個郵政編碼的索引，而不是郵政編碼本身。

來源

2012-01-08 CodeMed

我不認爲這是太辛苦與R.編寫爲R可接受的解決辦法？ – 2012-01-08 10:12:33

@RomanLuštrik，R可以工作。我希望看到示例代碼。 +1試圖幫助我。 – CodeMed 2012-01-08 20:13:36

您在ZipCodes上覆蓋myCSV，而不是'pts'SpatialPoints對象！錯誤消息是說「我不知道如何將數據幀覆蓋到SpatialPolygons數據框上」。覆蓋（pts，ZipCodes）應該可以工作。 – Spacedman 2012-01-09 22:44:27

Java庫做點/聚操作是JTS：

http://www.vividsolutions.com/jts/JTSHome.htm

你可能需要別的東西來閱讀shape文件，也許這樣的：

http://openmap.bbn.com/doc/api/com/bbn/openmap/layer/shape/ShapeFile.html

或可能java與OGR和GDAL的結合：

http://gdal.org/java/

然而，你可能只需要一個開源的GIS軟件包就可以做到這一點：Quantum GIS是我的最愛，但如果你想要一個基於Java的開源GIS，那麼OpenJump，gvSIG（www.osgeo.org對於這些東西的一切）。

在R，如果您已經閱讀了點座標爲矩陣或數據幀：

> xy 
      [,1]  [,2] 
[1,] 15.224275 -0.840066 
[2,] -1.207046 7.959644 
[3,] 9.397658 17.426323 
[4,] 28.242840 -29.581008 
[5,] 18.664603 15.361146 
[6,] 0.975846 8.534910 
[7,] -10.941825 10.438541 
[8,] -10.331097 20.260005 
[9,] 8.105570 9.595169 
[10,] -14.468177 14.366980

然後，使用maptools包及其依賴：

> require(maptools)

首先做出SpatialPoints對象從你的座標：

> pts = SpatialPoints(xy)

然後在你的shapefile中讀取：

> africa=readShapeSpatial("/path/to/africa.shp")

現在覆蓋：

> overlay(pts,africa) 
[1] 12 20 39 27 10 55 22 33 40 45

這是行號在shape文件的向量。您可以查找在shape文件很容易的屬性：

> africa$CNTRY_NAME[overlay(pts,africa)] 
[1] Congo  Ghana  Niger  Lesotho Chad  Togo  
[7] Guinea  Mauritania Nigeria Senegal 
61 Levels: Algeria Angola Benin Botswana Burkina Faso Burundi ... Zimbabwe

然後，如果你想要寫一個載體，一個CSV文件，

> write.csv(africa$CNTRY_NAME[overlay(pts,africa)],file="out.csv")

生產：

"","x" 
"1","Congo" 
"2","Ghana" 
"3","Niger" 
"4","Lesotho" 
"5","Chad" 
"6","Togo" 
"7","Guinea" 
"8","Mauritania" 
"9","Nigeria" 
"10","Senegal"

逗號分隔，引述，標題爲'x'。這些東西可以通過write.csv的其他選項進行調整。

如果您的點落在多邊形之外，返回覆蓋載體將是一個「NA」的值，你可能要測試這一點：

> if(any(is.na(overlay(pts,africa)))){stop("splash!")} 
Error: splash!

這份厚禮？

來源

2012-01-08 11:44:00 Spacedman

謝謝。 +1幫助我解決這個問題。我以前用R編碼過。如果您認爲這會使編碼變得簡單，我可以再次下載並研究它。我知道如何將.csv數據導入到R中，並且我可以很容易地將您的'匹配'輸出到.csv。但是您是否願意發佈您認爲我的R腳本需要完成的任務的更完整版本才能完成我在上面原始發佈中列出的整個任務？謝謝。 – CodeMed 2012-01-08 20:12:14

已編輯成完整的解決方案。 – Spacedman 2012-01-09 08:37:31

謝謝。我目前正在下載試圖嘗試這些方法。 +1幫助。 – CodeMed 2012-01-09 20:13:39