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

Divya
źródło
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
Możliwy duplikat - gis.stackexchange.com/questions/133099/…
klewis
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

kttii
źródło
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

Xunilk
źródło
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
T_Bacon
źródło