2013-08-28 61 views
0

我有一個csv文件,包含17,305個池塘的池塘面積和緯度和經度座標。對於每個池塘,我想確定1公里內所有池塘的座標。我是R新手,所以我想我可以適應一些最近的鄰居代碼。我發現這個環中的R預訂截止克拉雷:最近鄰居R代碼的適應性,以確定每個池塘1公里內的池塘位置

x<-runif(100) 
y<-runif(100) 

par(pty="s") 
plot(x,y,pch=16) 

distance<-function(x1, y1, x2, y2) sqrt((x2 − x1)^2 + (y2 − y1)^2) 

r<-numeric(100) 
nn<-numeric(100) 
d<-numeric(100) 
for (i in 1:100) { 
for (k in 1:100) d[k]<-distance(x[i],y[i],x[k],y[k]) 
r[i]<-min(d[-i]) 
nn[i]<-which(d==min(d[-i])) 
} 

for (i in 1:100) lines(c(x[i],x[nn[i]]),c(y[i],y[nn[i]])) 

我適應它和使用化石deg.dist函數使用所述半正矢式,而不是使用畢達哥拉斯的。

install.packages("fossil") 
library(fossil) 

Pond_A<-read.csv("C:\\ PondArea_data\\Pond_areas.csv") 

r<-numeric(17305) 
nn<-numeric(17305) 
d<-numeric(17305) 
for (i in 1:17305){ 
for (k in 1:17305) d[k]<-with(Pond_A,deg.dist(Longitude[i],Latitude[i],Longitude[k],Latitude[k])) 
    r[i]<-min(d[-i]) 
    nn<-which(d<=1) 
} 

這似乎給我所有的池塘在1公里的最後一個池塘的身份。但嘗試一下,因爲我可能無法弄清楚如何爲所有的池塘找到答案。如果有人能給我一個解決方案,並且可能解釋它爲什麼有效,我將非常感激。

感謝,

艾丹

+1

你看過'sp'包嗎?其中的'spDists'功能應該給你一些容易處理的東西。 – Frank

回答

0

您可以創建在rgeos包使用gWithinDistance布爾矩陣。 row/col值表示sp對象的rownames。然後,您可以將矩陣強制爲一個數據框並將其分配回sp對象。對於這個例子,我使用sp包中的meuse數據。

require(sp) 
require(rgeos) 
data(meuse) 
    coordinates(meuse) <- ~x+y 

# Create boolean matrix where TRUE is distance condition is |nnd <= d| TRUE else FALSE 
d=200 
DistMat <- gWithinDistance(meuse, meuse, dist=d, byid=TRUE) 

# Turn self-evaluation values to NA 
diag(DistMat) <- NA 

# Join back to data 
cids <- colnames(DistMat) 
    DistMat <- as.data.frame(DistMat) 
    names(DistMat) <- paste("NID", cids, sep=".") 
     [email protected] <- data.frame([email protected], DistMat) 
     str([email protected]) 
+0

非常感謝您的回答。對不起,它很長時間才能回覆我想確保我理解了答案並可以應用它。如果我使用緯度/長度座標,距離應該是度數,所以我需要花費更長的時間,因此我需要轉換爲東方和北方。 – user2358636