Wygładzanie wielokątów na mapie konturowej?

Oto mapa konturowa, dla której dostępne są wszystkie wielokąty poziomów.

Zapytajmy, jak wygładzić wielokąty, zachowując wszystkie wierzchołki w ich dokładnych lokalizacjach?

Rzeczywiście kontur jest nałożony na dane siatki, możesz więc zasugerować wygładzenie danych siatki, a tym samym powstały kontur będzie gładszy. Zauważ, że to nie działa zgodnie z moim pragnieniem, ponieważ funkcja wygładzania, taka jak filtr Gaussa, usunie małe paczki danych i zmieni zakres trzeciej zmiennej, np. Wysokość, która nie jest dozwolona w mojej aplikacji.

Właściwie szukam kawałka kodu (najlepiej w Pythonie ), który może wygładzać wielokąty 2D (dowolnego typu: wypukły, wklęsły, przecinający się itp.) W miarę bezbolesny (zapomnij o stronach kodów) i dokładny.

Do Twojej dyspozycji jest funkcja ArcGIS, która robi to doskonale, ale korzystanie z komercyjnych aplikacji innych firm nie jest moim wyborem.

1)

Scipy.interpolate:

Jak widać, powstałe splajny (czerwone) nie są zadowalające!

2)

Oto wynik przy użyciu kodu podanego tutaj . To nie działa dobrze!

3)

Dla mnie najlepszym rozwiązaniem powinno być coś takiego jak poniższy rysunek, w którym kwadrat jest wygładzany stopniowo poprzez zmianę tylko jednej wartości. Mam nadzieję na podobną koncepcję wygładzania dowolnej formy wielokątów.

Spełniając warunek, że splajn przechodzi przez punkty:

4)

Oto moja implementacja „pomysłu Whubera” wiersz po wierszu w Pythonie na jego danych. Możliwe są pewne błędy, ponieważ wyniki nie są dobre.

K = 2 to katastrofa, więc dla k> = 4.

5)

Usunąłem jeden punkt w problematycznej lokalizacji, a powstały splajn jest teraz identyczny z Whuberem. Ale wciąż pozostaje pytanie, dlaczego ta metoda nie działa we wszystkich przypadkach?

6)

Dobre wygładzanie danych Whubera może być następujące (rysowane przez oprogramowanie do grafiki wektorowej), w którym płynnie dodano dodatkowy punkt (porównaj z aktualizacją

4):

7)

Zobacz wynik z wersji Whubera w Pythonie dla niektórych ikonicznych kształtów:

Zauważ, że metoda wydaje się nie działać dla polilinii. Dla narożnej polilinii (konturu) chcę tego, co jest zielone, ale zrobiło się czerwone. Należy się tym zająć, ponieważ mapy konturowe są zawsze poliliniami, chociaż zamknięte polilinie można traktować jak wielokąty, jak w moich przykładach. Nie oznacza to również, że problem pojawiający się w aktualizacji 4 nie został jeszcze rozwiązany.

8) [mój ostatni]

Oto ostateczne rozwiązanie (nie idealne!):

Pamiętaj, że będziesz musiał coś zrobić z obszarem wskazanym przez gwiazdy. Być może w moim kodzie występuje błąd lub proponowana metoda wymaga dalszego rozwoju, aby uwzględnić wszystkie sytuacje i zapewnić pożądane wyniki.

Większość metod splajnowania sekwencji liczb splajnie wielokąty. Sztuczka polega na tym, aby splajny „zamykały się” płynnie w punktach końcowych. Aby to zrobić, owiń wierzchołki wokół końców. Następnie oddziel osobno współrzędne x i y.

Oto działający przykład w R. Wykorzystuje domyślną splineprocedurę sześcienną dostępną w pakiecie podstawowych statystyk. Aby uzyskać większą kontrolę, substytutem prawie każda procedura wolisz: Wystarczy upewnić się, dłutowanie it poprzez numery (czyli interpoluje je), a nie tylko przy użyciu ich jako „punkty kontrolne”.

#
# Splining a polygon.
#
#   The rows of 'xy' give coordinates of the boundary vertices, in order.
#   'vertices' is the number of spline vertices to create.
#              (Not all are used: some are clipped from the ends.)
#   'k' is the number of points to wrap around the ends to obtain
#       a smooth periodic spline.
#
#   Returns an array of points. 
# 
spline.poly <- function(xy, vertices, k=3, ...) {
    # Assert: xy is an n by 2 matrix with n >= k.

    # Wrap k vertices around each end.
    n <- dim(xy)[1]
    if (k >= 1) {
        data <- rbind(xy[(n-k+1):n,], xy, xy[1:k, ])
    } else {
        data <- xy
    }

    # Spline the x and y coordinates.
    data.spline <- spline(1:(n+2*k), data[,1], n=vertices, ...)
    x <- data.spline$x
    x1 <- data.spline$y
    x2 <- spline(1:(n+2*k), data[,2], n=vertices, ...)$y

    # Retain only the middle part.
    cbind(x1, x2)[k < x & x <= n+k, ]
}

Aby zilustrować jego użycie, stwórzmy mały (ale skomplikowany) wielokąt.

#
# Example polygon, randomly generated.
#
set.seed(17)
n.vertices <- 10
theta <- (runif(n.vertices) + 1:n.vertices - 1) * 2 * pi / n.vertices
r <- rgamma(n.vertices, shape=3)
xy <- cbind(cos(theta) * r, sin(theta) * r)

Spline to przy użyciu poprzedniego kodu. Aby wygładzić splajn, zwiększ liczbę wierzchołków od 100; aby było mniej gładkie, zmniejsz liczbę wierzchołków.

s <- spline.poly(xy, 100, k=3)

Aby zobaczyć wyniki, wykreślamy (a) pierwotny wielokąt w przerywanej czerwieni, pokazując odstęp między pierwszym a ostatnim wierzchołkiem (tj. Nie zamykając jego polilinii granicznej); oraz (b) splajn na szaro, jeszcze raz pokazując swoją szczelinę. (Ponieważ przerwa jest tak mała, jej punkty końcowe są podświetlone niebieskimi kropkami).

plot(s, type="l", lwd=2, col="Gray")
lines(xy, col="Red", lty=2, lwd=2)
points(xy, col="Red", pch=19)
points(s, cex=0.8)
points(s[c(1,dim(s)[1]),], col="Blue", pch=19)

Wiem, że to stary post, ale pojawił się w Google w poszukiwaniu czegoś, czego szukałem, więc pomyślałem, że opublikuję moje rozwiązanie.

Nie uważam tego za ćwiczenie dopasowania krzywej 2D, ale raczej ćwiczenie 3D. Uwzględniając dane jako 3D, możemy zapewnić, że krzywe nigdy się nie przecinają, i możemy wykorzystać informacje z innych konturów, aby poprawić nasze oszacowanie dla bieżącego.

Poniższy wyciąg iPython wykorzystuje interpolację sześcienną dostarczoną przez SciPy. Zauważ, że wartości Z, które narysowałem, nie są ważne, o ile wszystkie kontury mają jednakową odległość.

In [1]: %pylab inline
        pylab.rcParams['figure.figsize'] = (10, 10)
        Populating the interactive namespace from numpy and matplotlib

In [2]: import scipy.interpolate as si

        xs = np.array([0.0, 0.0, 4.5, 4.5,
                       0.3, 1.5, 2.3, 3.8, 3.7, 2.3,
                       1.5, 2.2, 2.8, 2.2,
                       2.1, 2.2, 2.3])
        ys = np.array([0.0, 3.0, 3.0, 0.0,
                       1.1, 2.3, 2.5, 2.3, 1.1, 0.5,
                       1.1, 2.1, 1.1, 0.8,
                       1.1, 1.3, 1.1])
        zs = np.array([0,   0,   0,   0,
                       1,   1,   1,   1,   1,   1,
                       2,   2,   2,   2,
                       3,   3,   3])
        pts = np.array([xs, ys]).transpose()

        # set up a grid for us to resample onto
        nx, ny = (100, 100)
        xrange = np.linspace(np.min(xs[zs!=0])-0.1, np.max(xs[zs!=0])+0.1, nx)
        yrange = np.linspace(np.min(ys[zs!=0])-0.1, np.max(ys[zs!=0])+0.1, ny)
        xv, yv = np.meshgrid(xrange, yrange)
        ptv = np.array([xv, yv]).transpose()

        # interpolate over the grid
        out = si.griddata(pts, zs, ptv, method='cubic').transpose()

        def close(vals):
            return np.concatenate((vals, [vals[0]]))

        # plot the results
        levels = [1, 2, 3]
        plt.plot(close(xs[zs==1]), close(ys[zs==1]))
        plt.plot(close(xs[zs==2]), close(ys[zs==2]))
        plt.plot(close(xs[zs==3]), close(ys[zs==3]))
        plt.contour(xrange, yrange, out, levels)
        plt.show()

Wyniki tutaj nie wyglądają najlepiej, ale przy tak małej liczbie punktów kontrolnych nadal są one całkowicie poprawne. Zwróć uwagę, jak wyciągana jest zielona dopasowana linia, aby podążać za szerszym niebieskim konturem.

Napisałem prawie dokładnie ten pakiet, którego szukasz ... ale był w Perlu i był ponad dekadę temu: GD :: Polyline . Używał 2D sześciennych krzywych Beziera i „wygładzał” dowolny wielokąt lub „polilinię” (wtedy nazywam to tak zwanym „LineString”).

Algorytm składał się z dwóch etapów: biorąc pod uwagę punkty w wielokącie, dodaj dwa punkty kontrolne Beziera między każdym punktem; następnie wywołaj prosty algorytm, aby dokonać częściowego przybliżenia splajnu.

Druga część jest łatwa; pierwsza część była trochę sztuką. Tutaj był wgląd: rozważyć „Segment kontrolny” wierzchołek n: vN. Segment kontrolny był trzy współliniowe punkty: [cNa, vN, cNb]. Punktem centralnym był wierzchołek. Nachylenie tego kontrolnego odcinka było równe nachyleniu od wierzchołka N-1 do wierzchołka N + 1. Długość lewej części tego odcinka wynosiła 1/3 długości od wierzchołka N-1 do wierzchołka N, a długość prawej części tego odcinka wynosiła 1/3 długości od wierzchołka N do wierzchołka N + 1.

Jeśli oryginalna krzywa była cztery wierzchołki: [v1, v2, v3, v4]wtedy każdy wierzchołek teraz dostać segment kontrolny w postaci: [c2a, v2, c2b]. Złóż je razem w ten sposób: [v1, c1b, c2a, v2, c2b, c3a, v3, c3b, c4a, v4]i przełam je po cztery jako cztery punkty Beziera: [v1, c1b, c2a, v2]a następnie [v2, c2b, c3a, v3]i tak dalej. Ponieważ [c2a, v2, c2b]były współliniowe, wynikowa krzywa będzie gładka w każdym wierzchołku.

Spełnia to zatem również wymagania dotyczące parametryzacji „szczelności” krzywej: użyj mniejszej wartości niż 1/3 dla „ciasniejszej” krzywej, większej dla dopasowania „zapętlonego”. W obu przypadkach wynikowa krzywa zawsze przechodzi przez pierwotne podane punkty.

Spowodowało to gładką krzywą „opisującą” oryginalny wielokąt. Miałem też sposób na „wpisanie” gładkiej krzywej ... ale nie widzę tego w kodzie CPAN.

W każdym razie nie mam obecnie dostępnej wersji w Pythonie ani żadnych liczb. ALE ... jeśli / kiedy przesyłam to do Pythona, z pewnością opublikuję tutaj.