2012-03-20 13 views
10

Próbuję znaleźć odległość ortogonalną między zbiorem współrzędnych lokalizacji a zbiorem linii (dróg lub rzek). Zestaw punktów ma postać par szerokości/długości geograficznej, a linie są w pliku kształtu (.shp). Umieszczenie ich na mapie nie stanowi problemu, używając albo maptools lub PBSmapping. Ale moim podstawowym problemem jest znalezienie minimalnej odległości, którą trzeba pokonać, aby dotrzeć do drogi lub rzeki. Czy jest jakiś sposób to zrobić w R?R/GIS: Znajdź ortogonalny dystans między położeniem a najbliższą linią

Odpowiedz

19

Jeśli dobrze rozumiem, możesz to zrobić po prostu z gDistance w pakieciew pakiecie.

przeczytać w linii jak i punktów SpatialLines/DataFrame jak SpatialPoints/DataFrame a następnie w pętli nad każdym punktem obliczaniu odległości za każdym razem:

require(rgeos) 
## untested code 
shortest.dists <- numeric(nrow(sp.pts)) 
for (i in seq_len(nrow(sp.pts)) { 
    shortest.dists[i] <- gDistance(sp.pts[i,], sp.lns) 
} 

Tutaj sp.pts to punkty obiektu przestrzennego, a sp.lns to linie obiektu przestrzennego.

Musisz pętli, aby porównać tylko jedną współrzędną w sp.pts z całością wszystkich geometrii linii w sp.lns, w przeciwnym razie otrzymasz odległość od wartości zagregowanej we wszystkich punktach.

Ponieważ twoje dane znajdują się w szerokości/długości geograficznej, powinieneś przekształcić obie linie i punkty w odpowiednią projekcję, ponieważ funkcja gDistance przyjmuje odległość kartezjańską.

więcej dyskusji i przykład (edit)

Byłoby miłe dostać najbliższego punktu na linii/s, a nie tylko od odległości, ale to otwiera inną opcję, która jest czy trzeba najbliższy koordynować wzdłuż linii lub rzeczywistego skrzyżowania z odcinkiem linii, który jest bliższy niż jakikolwiek istniejący wierzchołek. Jeśli twoje wierzchołki są wystarczająco gęste, aby różnica nie ma znaczenia, użyj pakietu spDistsN1 w pakiecie . Musisz wyodrębnić wszystkie współrzędne z każdej linii w zestawie (nie twarde, ale trochę brzydkie), a następnie przelotować nad każdym punktem zainteresowania, obliczając odległość do wierzchołków linii - wtedy możesz znaleźć, który jest najkrótszy i wybrać tę współrzędną z zestawu wierzchołków, dzięki czemu można łatwo uzyskać odległość i współrzędne. Nie trzeba również projektować, ponieważ funkcja może używać odległości elipsoidalnych z argumentem longlat = TRUE.

library(maptools) 

## simple global data set, which we coerce to Lines 
data(wrld_simpl) 

wrld_lines <- as(wrld_simpl, "SpatialLinesDataFrame") 

## get every coordinate as a simple matrix (scary but quick) 
wrld_coords <- do.call("rbind", lapply([email protected], function(x1) do.call("rbind", lapply([email protected], function(x2) [email protected][-nrow([email protected]), ])))) 

Sprawdź to interaktywnie, musisz to zmienić, aby zapisać współrzędne lub minimalne odległości. Spowoduje to wydrukowanie linii i poczekanie, aż klikniesz w dowolnym miejscu wykresu, a następnie narysuje linię od twojego kliknięcia do najbliższego wierzchołka wierzchołka na linii.

## no out of bounds clicking . . . 
par(mar = c(0, 0, 0, 0), xaxs = "i", yaxs = "i") 

plot(wrld_lines, asp = "") 

n <- 5 

for (i in seq_len(n)) { 
xy <- matrix(unlist(locator(1)), ncol = 2) 
    all.dists <- spDistsN1(wrld_coords, xy, longlat = TRUE) 
    min.index <- which.min(all.dists) 
    points(xy, pch = "X") 
lines(rbind(xy, wrld_coords[min.index, , drop = FALSE]), col = "green", lwd = 2) 
} 
+0

Nicea i pouczające odpowiedź. Dziękuję za to. Wygląda na to, że 'rgeos' nie ma żadnej funkcji, która zwraca (lub) najbliższą lokalizację. Czy zdarzyło ci się mieć sugestię dotyczącą funkcji, która to robi?Wygląda na to, że 'project2segment()' w 'spatstat' to robi, ale niestety nie używa obiektów klasy' sp', które byłyby poręczniejsze ... –

+0

Będę aktualizował, myślałem o tym więcej i są kilka subtelności. – mdsumner

+0

Konwersja między spatstat i sp nie jest trudna, istnieją narzędzia wymuszające w maptools, takie jak spst <- as.psp (wrld_lines) – mdsumner

3

Pakiet geosphere posiada funkcję dist2line że robi to dla lon/danych lat. Może używać obiektów lub macierzy przestrzennych *.

line <- rbind(c(-180,-20), c(-150,-10), c(-140,55), c(10, 0), c(-140,-60)) 
pnts <- rbind(c(-170,0), c(-75,0), c(-70,-10), c(-80,20), c(-100,-50), 
     c(-100,-60), c(-100,-40), c(-100,-20), c(-100,-10), c(-100,0)) 
d <- dist2Line(pnts, line) 
d 

Ilustracja wyników

plot(makeLine(line), type='l') 
points(line) 
points(pnts, col='blue', pch=20) 
points(d[,2], d[,3], col='red', pch='x') 
for (i in 1:nrow(d)) lines(gcIntermediate(pnts[i,], d[i,2:3], 10), lwd=2) 
Powiązane problemy