Algoritmo para averiguar los puntos de inflexión de una polilínea

Question

Intento averiguar los puntos de inflexión, es decir, los puntos en los que empiezan y terminan las curvas de una línea. Si te fijas en la imagen, la línea verde puede ser una carretera o un arroyo, y los puntos negros son los puntos donde empiezan y acaban las curvas.

¿Cuáles serían los pasos de alto nivel para automatizar la generación de estos puntos? Tengo ArcGIS desktop y soy bastante hábil con ArcObjects.

Answer 1

Cuando la curva está formada por segmentos de línea, todos los puntos interiores de esos segmentos son puntos de inflexión, lo que no es interesante. En su lugar, debe pensarse que la curva se aproxima por la vértices de esos segmentos. Si se empalma una curva a trozos dos veces diferenciable a través de esos segmentos, se puede calcular la curvatura. Un punto de inflexión, estrictamente hablando, es un lugar donde la curvatura es cero.

En el ejemplo hay tramos largos en los que la curvatura es casi nula. Esto sugiere que los puntos indicados deberían aproximarse a los extremos de esos tramos de regiones de baja curvatura.

Por lo tanto, un algoritmo eficaz dividirá los vértices, calculará la curvatura a lo largo de un conjunto denso de puntos intermedios, identificará los rangos de curvatura cercana a cero (utilizando alguna estimación razonable de lo que significa estar "cerca") y marcará los puntos finales de esos rangos.

Aquí está trabajando R para ilustrar estas ideas. Empecemos con una cadena de líneas expresada como una secuencia de coordenadas:

xy <- matrix(c(5,20, 3,18, 2,19, 1.5,16, 5.5,9, 4.5,8, 3.5,12, 2.5,11, 3.5,3, 
               2,3, 2,6, 0,6, 2.5,-4, 4,-5, 6.5,-2, 7.5,-2.5, 7.7,-3.5, 6.5,-8), ncol=2, byrow=TRUE)

Spline el x y y coordenadas por separado para conseguir una parametrización de la curva. (El parámetro se llamará time .)

n <- dim(xy)[1]
fx <- splinefun(1:n, xy[,1], method="natural")
fy <- splinefun(1:n, xy[,2], method="natural")

Interpolar las splines para el trazado y el cálculo:

time <- seq(1,n,length.out=511)
uv <- sapply(time, function(t) c(fx(t), fy(t)))

Necesitamos una función para calcular la curvatura de una curva parametrizada. Necesita estimar las derivadas primera y segunda del spline. Con muchos splines (como los cúbicos) esto es un cálculo algebraico fácil. R proporciona las tres primeras derivadas automáticamente. (En otros entornos, es posible que se quiera calcular las derivadas numéricamente).

curvature <- function(t, fx, fy) {
  # t is an argument to spline functions fx and fy.
  xp <- fx(t,1); yp <- fy(t,1)            # First derivatives
  xpp <- fx(t,2); ypp <- fy(t,2)          # Second derivatives
  v <- sqrt(xp^2 + yp^2)                  # Speed
  (xp*ypp - yp*xpp) / v^3                 # (Signed) curvature
  # (Left turns have positive curvature; right turns, negative.)
}

kappa <- abs(curvature(time, fx, fy))     # Absolute curvature of the data

Propongo estimar un umbral de curvatura cero en cuanto a la extensión de la curva. Este es, al menos, un buen punto de partida; debería ajustarse en función de la tortuosidad de la curva (es decir, aumentar para las curvas más largas). Esto se utilizará más tarde para colorear los gráficos según la curvatura.

curvature.zero <- 2*pi / max(range(xy[,1]), range(xy[,2])) # A small threshold
i.col <- 1 + floor(127 * curvature.zero/(curvature.zero + kappa)) 
palette(terrain.colors(max(i.col)))                        # Colors

Ahora que los vértices han sido empalmados y la curvatura calculada, sólo queda encontrar los puntos de inflexión . Para mostrarlos podemos trazar los vértices, trazar la spline y marcar los puntos de inflexión en ella.

plot(xy, asp=1, xlab="x",ylab="y", type="n")
tmp <- sapply(2:length(kappa), function(i) lines(rbind(uv[,i-1],uv[,i]), lwd=2, col=i.col[i]))
points(t(sapply(time[diff(kappa < curvature.zero/2) != 0], 
       function(t) c(fx(t), fy(t)))), pch=19, col="Black")
points(xy)

Los puntos abiertos son los vértices originales en xy y los puntos negros son los puntos de inflexión identificados automáticamente con este algoritmo. Dado que la curvatura no puede calcularse de forma fiable en los puntos finales de la curva, esos puntos no están especialmente marcados.

Answer 2

Puede utilizar el Herramienta de densificación . Para este caso, se elige densificar por ángulo, A continuación, elija el ángulo máximo aceptado en una línea recta. A continuación, aplique a la línea resultante a la herramienta Dividir la línea en los vértices . Por último, elimine las líneas que tengan shape_length menor que la longitud mínima de la carretera.

En esta imagen, vemos tres pasos:

1- Densificar la línea utilizando el ángulo. He utilizado 10 grados como parámetro, y hemos utilizado splitline. En la imagen, la línea curva está en su fase inicial.

arcpy.Densify_edit("line" , "ANGLE" , "","",10)
arcpy.SplitLine_management("line" , "line_split")

2- Seleccionar los segmentos que tengan shape_length no redundante. Como vemos en la tabla, no he seleccionado esas longitudes redundantes. Entonces, los selecciono en una nueva clase de característica.

arcpy.Select_analysis("line_split" , "line_split_selected")

3- Hemos extraído los vértices situados en las aristas de las líneas, que son puntos de inflexión.

arcpy.FeatureVerticesToPoints_management("line_split_selected" , "line_split_pnt" , "DANGLE")

Answer 3

Puede utilizar el Generalizar que tiene como parámetro el desplazamiento máximo respecto a la línea original, por lo que puede elegir el desplazamiento que se ajuste a su caso.

Si llamamos a la línea original "line_cur", y a la generalizada "line_gen", entonces podríamos recortar "line_cur" por "line_gen". El resultado será el segmento recto de "line_cur". Luego podríamos limpiar algunos segmentos muy cortos borrándolos con una consulta sql que seleccione el Shape_length mayor que la longitud mínima del camino.