Algorytm znajdowania punktów przegięcia dla polilinii

Próbuję znaleźć punkty przegięcia, tj. Punkty, w których krzywe linii zaczynają się i kończą. Jeśli spojrzysz na obraz, zielona linia może być drogą lub strumieniem, a czarne punkty to punkty, w których krzywe zaczynają się i kończą.

Jakie byłyby kroki na wysokim szczeblu w celu automatyzacji generowania tych punktów? Mam pulpit ArcGIS i jestem całkiem przydatny w ArcObjects.

Gdy krzywa składa się z odcinków linii, wówczas wszystkie punkty wewnętrzne tych odcinków są punktami przegięcia, co nie jest interesujące. Zamiast tego należy uważać, że krzywa jest zbliżona do wierzchołków tych segmentów. Poprzez rozcięcie częściowo segmentowanej krzywej o dwóch różnicach między tymi segmentami, możemy następnie obliczyć krzywiznę. Mówiąc wprost, punkt przegięcia to miejsce, w którym krzywizna wynosi zero.

W tym przykładzie występują długie odcinki, w których krzywizna wynosi prawie zero. Sugeruje to, że wskazane punkty powinny przybliżać końce takich odcinków regionów o niskiej krzywiźnie.

Dlatego skuteczny algorytm splajtuje wierzchołki, oblicza krzywiznę wzdłuż gęstego zestawu punktów pośrednich, identyfikuje zakresy prawie zerowej krzywizny (wykorzystując rozsądne oszacowanie, co to znaczy być „blisko”), i zaznacza punkty końcowe tych zakresów .

Oto działający Rkod ilustrujący te pomysły. Zacznijmy od ciągu linii wyrażonego jako ciąg współrzędnych:

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 x i y współrzędne osobno do uzyskania parametryzacji krzywej. (Parametr zostanie wywołany time.)

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

Interpoluj splajny na potrzeby kreślenia i obliczeń:

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

Potrzebujemy funkcji do obliczenia krzywizny sparametryzowanej krzywej. Musi oszacować pierwszą i drugą pochodną splajnu. W przypadku wielu splajnów (takich jak splajny sześcienne) jest to łatwe obliczenie algebraiczne. Rzapewnia trzy pierwsze pochodne automatycznie. (W innych środowiskach można chcieć obliczać pochodne numerycznie).

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

Proponuję oszacować próg zerowej krzywizny pod względem zasięgu krzywej. To przynajmniej dobry punkt wyjścia; należy go wyregulować zgodnie z krętością krzywej (to znaczy zwiększyć dla dłuższych krzywych). Będzie to później wykorzystane do pokolorowania wykresów zgodnie z krzywizną.

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

Teraz, gdy wierzchołki zostały spline i obliczono krzywiznę, pozostaje tylko znaleźć punkty przegięcia . Aby je pokazać, możemy narysować wierzchołki, narysować splajn i zaznaczyć na nim punkty przegięcia.

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)

Punkty otwarte to oryginalne wierzchołki, xya czarne punkty to punkty przegięcia automatycznie identyfikowane za pomocą tego algorytmu. Ponieważ krzywizny nie można wiarygodnie obliczyć w punktach końcowych krzywej, punkty te nie są specjalnie oznaczone.

Możesz użyć narzędzia Densify . W tym przypadku wybierz zagęszczenie według kąta, następnie wybierz maksymalny kąt akceptowany w linii prostej. Następnie zastosuj do linii wyników do narzędzia Linia podziału w wierzchołkach . Na koniec usuń linie o kształcie_długości mniejszym niż minimalna długość drogi.

Na tym zdjęciu widzimy trzy kroki:

1- Zagęścić linię za pomocą kąta. Jako parametru użyłem 10 stopni i zastosowaliśmy linię podziału. Na zdjęciu zakrzywiona linia znajduje się w początkowej fazie.

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

2 - Wybierz segmenty, w których długość kształtu nie jest redundantna. Jak widać w tabeli, nie wybrałem tych zbędnych długości. Następnie wybieram je do nowej klasy funkcji.

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

3- Wyodrębniliśmy wierzchołki znajdujące się na krawędziach linii, które są punktami przegięcia.

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

Możesz użyć narzędzia Uogólnij , który ma maksymalne przesunięcie od oryginalnej linii jako parametr, dzięki czemu możesz wybrać przesunięcie, które pasuje do twojego przypadku.

Jeśli nazwiemy pierwotną linię „line_cur”, a uogólnioną „line_gen”, możemy przyciąć „line_cur” przez „line_gen”. Wynikiem będzie prosty segment „line_cur”. Następnie możemy wyczyścić bardzo krótki odcinek, usuwając go za pomocą zapytania SQL, które wybiera Shape_length większą niż minimalna długość drogi.