Odległość między środkiem ciężkości a najdalszym punktem wielokąta


12

Mam wiejską warstwę wielokąta, która ma ponad 6 000 000 rekordów. Obliczyłem centroid każdej wioski. Chcę znaleźć odległość między środkiem ciężkości a najdalszym węzłem każdego wielokąta. Sprawdź zdjęcie poniżej w celach informacyjnych. Czarne linie to granice wielokątów. wprowadź opis zdjęcia tutaj

wprowadź opis zdjęcia tutaj

wprowadź opis zdjęcia tutaj


ciekawe ... Właśnie zrobiłem ten piątek z postgis, aby stworzyć okrąg wokół wielokąta. Potrzebuję kilku minut, aby poszukać kodu, którego użyłem .. i.stack.imgur.com/EKnkg.png
— kttii

1
Najpierw możemy potrzebować wiedzieć, jakie programy masz do dyspozycji. Jak również stworzyłeś te centroidy i węzły? (Nawet jeśli wydaje się to nieco oczywiste, że węzły na wielokątach służą do wyznaczania granic twoich kształtów, ale czy dodałeś dodatkowy punkt na nich?)
— Moreau Colin

Czy lokalizacja środka ciężkości jest ważna? Jak je stworzyłeś?
— GISGe


Jeśli środek ciężkości jest naprawdę centralny, to promień najmniejszego koła wyśrodkowanego w tym punkcie pasuje do wielokąta ( en.wikipedia.org/wiki/Smallest-circle_problem )
— Mark Ireland

Odpowiedzi:



15

Korzystając z PostGIS, użyłem ST_ConvexHull, aby uprościć wielokąt i uzyskać szybszy wynik:

Zdobądź najdalszy punkt:

SELECT Villages_v4_Trial_region.geom as FarPoint from (
SELECT ST_PointN(ST_ExteriorRing(ST_ConvexHull(Villages_v4_Trial_region.geom)),
generate_series(1, ST_NPoints(ST_ExteriorRing(ST_ConvexHull(Villages_v4_Trial_region.geom))))) as points, 
geom
FROM Villages_v4_Trial_region
ORDER BY ST_MaxDistance(points,ST_Centroid(Villages_v4_Trial_region.geom)) DESC
LIMIT 1;

A jeśli chcesz utworzyć okrąg z centrum ciężkości:

SELECT ST_Buffer(Center,ST_Distance(Center,FarPoint)) as Circle
FROM (
SELECT Villages_v4_Trial_region.geom as FarPoint, Center from (
    SELECT ST_PointN(ST_ExteriorRing(ST_ConvexHull(Villages_v4_Trial_region.geom)),
    generate_series(1, ST_NPoints(ST_ExteriorRing(ST_ConvexHull(Villages_v4_Trial_region.geom))))) as points,
    ST_Centroid(Villages_v4_Trial_region.geom) as Center, 
    geom
    FROM Villages_v4_Trial_region
    ) as Villages_v4_Trial_region
    ORDER BY ST_MaxDistance(points,Center) DESC
    LIMIT 1) as foo;

wprowadź opis zdjęcia tutaj


Prosty, szybki, wydajny. Dziękujemy za opublikowanie tego, ponieważ pomoże mi to również w tym, nad czym teraz pracuję.
— Moreau Colin

@kttii Nie wiem, jak korzystać z PostGIS. Czy możesz podać prostsze rozwiązanie, w trybie łukowym, mapinfo lub qgis
— Divya

@kttii Więc zainstalowałem Postgresql. Skopiowałem dokładnie to zapytanie, ale dało BŁĄD: kolumna „the_geom” nie istnieje. Co ja robię?
— Divya

the_geom należy zastąpić nazwą pola geometrii. Będziesz musiał także umieścić swoje dane w PostgreSQL. PostgreSQL to baza danych podobna do MSSQL. PostGIS jest rozszerzeniem, dzięki któremu baza danych jest przestrzennie świadoma i zapewnia wszystkie funkcje ST_.
— kttii

@kttii Zaktualizowałem nazwę pola z the_geom do „gid” w mojej bazie danych. Po ponownym uruchomieniu zapytania otrzymałem ten BŁĄD: funkcja st_convexhull (liczba całkowita) nie istnieje
— Divya

4

Za pomocą następnego kodu PyQGIS :

from math import sqrt

layer = iface.activeLayer()

feats = [ feat for feat in layer.getFeatures() ]

n = len(feats)

centroids = [ feat.geometry().centroid().asPoint() for feat in feats ]
polygons = [ feat.geometry().asPolygon()[0] for feat in feats ]

lengths = []

for i, pol in enumerate(polygons):
    max_dist = 0
    idx_j = 0
    for j, point in enumerate(pol):
        dist = sqrt(centroids[i].sqrDist(point))
        if dist > max_dist:
            max_dist = dist
            idx_j = j
    print i, idx_j, max_dist
    lengths.append([centroids[i], pol[idx_j]])

crs = layer.crs()
epsg = crs.postgisSrid()

uri = "LineString?crs=epsg:" + str(epsg) + "&field=id:integer""&index=yes"

mem_layer = QgsVectorLayer(uri,
                           'max_distance',
                           'memory')

prov = mem_layer.dataProvider()

feats = [ QgsFeature() for i in range(n) ]

for i, feat in enumerate(feats):
    feat.setAttributes([i])
    feat.setGeometry(QgsGeometry.fromPolyline(lengths[i]))

prov.addFeatures(feats)

QgsMapLayerRegistry.instance().addMapLayer(mem_layer)

i ten plik kształtu (z 11 funkcjami):

wprowadź opis zdjęcia tutaj

Mam warstwę pamięci, w której polilinie były odległościami między środkiem ciężkości a najdalszym punktem każdego wielokąta (cecha); jak widać na następnym obrazku:

wprowadź opis zdjęcia tutaj

W konsoli Pythona QGIS wydrukowano również indeks elementu, indeks punktu w elemencie, w którym odległość od środka ciężkości jest maksymalna, a na koniec maksymalna.

wprowadź opis zdjęcia tutaj


Nie wiem, jak korzystać z PyQGIS. Czy możesz podać prostsze rozwiązanie, w trybie łukowym, mapinfo lub qgis
— Divya

1
Wypróbuj ten link, aby uzyskać pomoc w rozpoczęciu korzystania z PyQgis spatialgalaxy.net/2014/10/09/…
— kttii

0

Wygląda na to, że korzystasz z MapInfo, oto funkcja MapBasic, którą napisałem jakiś czas temu dla wewnętrznego narzędzia, nad którym pracowałem. Jako argument przyjmuje węzeł źródłowy (punkt środkowy) i obiekt regionu (wielokąt) i zwraca obiekt punktowy w najdalszym węźle wielokąta od punktu źródłowego.

Function GetFurthest(ByVal oNode1 as Object, ByVal oObj as Object) as Object

Dim sourceE,sourceN,East,North,Longest,Dist as Float,
    nNodes,nPolys,i,j as SmallInt,
    oNode2 as Object

    sourceE = CentroidX(oNode1)
    sourceN = CentroidY(oNode1)
    Longest = 0

    nPolys = ObjectInfo(oObj,OBJ_INFO_NPOLYGONS)
    For i = 1 to nPolys
        nNodes = ObjectInfo(oObj,OBJ_INFO_NPOLYGONS+nPolys)
        For j = 1 to nNodes
            East = ObjectNodeX(oObj,i,j)
            North = ObjectNodeY(oObj,i,j)
            Dist = Distance(sourceE,sourceN,East,North,"m")
            If Dist > Longest then
                Longest = Dist
                oNode2 = CreatePoint(East,North)
            End if
        Next
    Next

    GetFurthest = oNode2

End Function
Korzystając z naszej strony potwierdzasz, że przeczytałeś(-aś) i rozumiesz nasze zasady używania plików cookie i zasady ochrony prywatności.
Licensed under cc by-sa 3.0 with attribution required.