R/GIS: Encuentre la distancia ortogonal entre una ubicación y la línea más cercana

Question

Estoy tratando de encontrar la distancia ortogonal entre un conjunto de coordenadas de ubicación y un conjunto de líneas (carreteras o ríos). El conjunto de puntos tiene forma de pares de latitud/longitud y las líneas están en un archivo de forma (.shp). Trazarlos en un mapa no es un problema, usando maptools o PBSmapping. Pero mi problema básico es encontrar la distancia mínima que uno tiene que recorrer desde una ubicación para llegar a una carretera o un río. ¿Hay alguna manera de hacer esto en R?

Answer 1

Si entiendo correctamente, puede hacer esto simplemente con gDistance en el paquete rgeos.

Leer en las líneas como SpatialLines/DataFrame y puntos como SpatialPoints/DataFrame y luego bucle sobre cada punto de cálculo de la distancia cada vez:

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) 
} 

Aquí sp.pts es objeto los puntos espaciales, y sp.lns es objeto de las líneas espacial.

Debe realizar un ciclo para que solo compare una sola coordenada en sp.pts con la totalidad de las geometrías de todas las líneas en sp.lns, de lo contrario, obtendrá la distancia desde un valor agregado en todos los puntos.

Dado que sus datos están en latitud/longitud, debe transformar tanto las líneas como los puntos en una proyección adecuada, ya que la función gDistance asume una distancia cartesiana.

MÁS DISCUSIÓN Y EJEMPLO (editar)

Sería ordenada para conseguir el punto de la línea más próxima/s en lugar de sólo la distancia, pero esto abre otra opción, que es si usted necesita el más cercana coordenadas a lo largo de una línea, o una intersección real con un segmento de línea que está más cerca que cualquier vértice existente. Si sus vértices son lo suficientemente densos como para que la diferencia no importe, entonces use spDistsN1 en el paquete sp. Tendría que extraer todas las coordenadas de cada línea del conjunto (no duro, pero un poco feo) y luego recorrer cada punto de interés calculando la distancia a los vértices de la línea, entonces puede encontrar cuál es el más corto y seleccionar esa coordenada del conjunto de vértices, para que pueda tener la distancia y la coordenada fácilmente. No es necesario proyectar tampoco, ya que la función puede usar distancias elipsoidales con el argumento 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(wrld_lines@lines, function(x1) do.call("rbind", lapply(x1@Lines, function(x2) x2@coords[-nrow(x2@coords), ])))) 

Compruébelo de forma interactiva, tendrá que modificar esto para guardar los coords o las distancias mínimas. Esto trazará las líneas y esperará a que haga clic en cualquier lugar de la gráfica, luego dibujará una línea desde su clic hasta el vértice más cercano en una línea.

## 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) 
} 

Answer 2

El paquete geosphere tiene la función dist2line que hace esto por lon/datos lat. Puede usar objetos o matrices Espaciales *.

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 

Ilustración de los resultados

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)