Wie man dataframe Zeilen, in denen X-und Y-Koordinaten außerhalb des Polygons

Ich versuche zu Adresse das folgende Problem. Nehmen wir einen dataframe (geladen aus einer txt-Datei) mit der folgenden Struktur (und Tausende von Zeilen):

foo.head()

         X            Y       Z 
 0  125417.5112  536361.8752 -1750.0
 1  127517.7647  533925.8644 -1750.0
 2  128144.1000  533199.4000 -1750.0
 3  128578.8385  532904.9288 -1750.0
 4  125417.5112  536361.8752 -1750.0
 ....

Die Daten darstellt, X -, Y-und Z-Koordinaten.

Ich auch eine Menge Punkte, die definieren, die ein geschlossenes polygon. Diese sind in ein numpy-array:

polypoints

array([[ 125417.5112,  536361.8752],
       [ 127517.7647,  533925.8644],
       [ 128144.1   ,  533199.4   ],
       ....
       [ 125417.5112,  536361.8752]])

Wie kann ich meine filter dataframe so löschen Sie die Zeilen, die NICHT fallen in dem geschlossenen polygon?

Ich habe versucht die Definition des Polygons mit shapely.geometry polygon. By doing:

poly = Polygon(polypoints)

Diese funktioniert einwandfrei. Aber ich bin an einem Verlust, wie Sie weiter mit diesem.

Hilfe ist sehr willkommen

– – – – – EDIT – – – – –
Bitte siehe unten für die aktualisierte Lösung

Der klassische Algorithmus für diese zu zeichnen ist eine Linie vom Punkt zur Unendlichkeit, der keine Schnittmenge mit irgendeinem der Punkte und zählen, wie viele Kanten kreuzt. Seltsam für drinnen, noch für draußen.
Haben Sie einen Blick auf Geopandas: geopandas.org/set_operations.html
das wäre für eine „einfache“ geometrische Form, denke ich. In meinem Fall ist es ein Komplexes polygon, damit es nicht so geradlinig. Auch ich bin auf der Suche nach einer Lösung, die arbeitet die ganze Zeit
die overlay-Funktion sieht wirklich vielversprechend aus. In Sie suchen im moment. Brauchen, um herauszufinden, ob es angewendet werden kann, auch mit polygon-und Punkte
Es funktioniert mit jedem Komplexität polygon, es wird einfach teurer zu berechnen als die Anzahl der Kanten geht.

InformationsquelleAutor Red Sparrow | 2018-02-09

Tags: dataframe, geospatial, pandas, point-in-polygon, python-3.x

3 Replies

2

Die original Lösung von @MrT funktioniert Super. Dennoch, mit Blick auf geopandas wie vorgeschlagen von @Rutger Kassies, ich habe auch eine andere Lösung gefunden. Zuerst braucht man zum installieren der geopandas Paket. Dann wird der folgende code funktioniert bei mir:
```
import geopandas as gpd
from shapely.geometry import Point, Polygon, MultiPolygon
# load the data that should be cropped by the polygon
# this assumes that the csv file already includes 
# a geometry column with point data as performed below
dat_gpd = gpd.GeoDataFrame.from_csv(r'{}\data_to_crop.csv'.format(savedir), sep='\t')

# load the data of the polygon as a dataframe
arr_df = pd.DataFrame(data, columns=['X','Y','Z'])

# make shapely points out of the X and Y coordinates
point_data = [Point(xy) for xy in zip(arr_df.X, arr_df.Y)]

# assign shapely points as geometry to a geodataframe
# Like this you can also inspect the individual points if needed
arr_gpd = gpd.GeoDataFrame(arr_df, geometry=point_data)

# define a shapely polygon from X and Y coordinates of the shapely points
polygo = Polygon([[p.x, p.y] for p in arr_gpd.geometry])

# assing defined polygon to a new dataframe
pol_gpd= gpd.GeoDataFrame()
pol_gpd['geometry'] = None
pol_gpd.loc[0,'geometry'] = polygo

# define a new dataframe from the spatial join of the dataframe with the data to be cropped
# and the dataframe with the polygon data, using the within function.
dat_fin = gpd.sjoin(dat_gpd, pol_gpd, op = 'within')
```
Hoffe, das hilft, wenn jemand vor einer ähnlichen problem. Auch, weitere Infos über die räumliche Verknüpfung kann gefunden werden auf der website geopandas. Beachten Sie, dass diese Funktionalität nicht benötigen eine operation zwischen Polygonen, funktioniert aber auch mit Punkte und Polygone

–EDIT —
```
%timeit gpd.sjoin(dat_gpd, pol_gpd, op = 'within')
31.8 s ± 108 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

%timeit dat_gpd['inpoly'] = dat_gpd.apply(lambda row: polygo.intersects(Point(row["X"], row["Y"])), axis = 1)
1min 26s ± 389 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)
```
Scheint es, dass die geo-pandas-Funktion ist viel schneller. Aber um fair zu sein die nicht-geo pandas Lösung hat auch zum konvertieren der X-und Y formschöne zeigen Sie Elemente und führen Sie dann die Kreuzung Bewertung
- Vielen Dank für das hinzufügen Ihrer eigenen Lösung. Ich hatte das Gefühl, dass formschöne und pandas sollten arbeiten gut zusammen, da beide basieren auf numpy. Ein wenig enttäuschend ist jedoch, dass Sie noch zu berechnen, die jedes polygon einzeln. Sie können die Geschwindigkeit für beide Ansätze mit der timeit – Modul für Daten-set von mittlerer Größe eine Vorstellung zu bekommen, welcher Weg Sie sparen Zeit, wenn Sie haben eine Menge Daten zu STEMMEN. +1
- Ich habe die timings. Es ist ein Unterschied, aber ich bin nicht sicher, ich bin den Vergleich äpfel mit äpfeln hier. Auch der code, den ich zur Verfügung gestellt werden, möglicherweise nicht der beste Weg, um zu definieren, formschöne Elemente in geopandas. Es funktioniert aber vielleicht nicht die eleganteste. Vielleicht kann mir jemand helfen, damit es mehr pythonic
- Es ist zu erwarten, dass eine echte geopandas code schneller arbeitet. Ich bin nicht überrascht. Für eine Leistungssteigerung, Sie könnten Ihren code unter die Lupe genommen und auf CodeReview.
InformationsquelleAutor Red Sparrow
2

Ich bin nicht so vertraut mit shapely. Vielleicht haben Sie eine echte pandas zu unterstützen. Afaik unterstützen Sie vektorisiert numpy-Funktionen, so wäre ich nicht überrascht.

Ein Weg, um herauszufinden, den die Punkte innerhalb eines gegebenen Polygons, wäre die Verwendung von pandas apply() Funktion:
```
import pandas as pd
from shapely.geometry import Polygon, Point
#your dataframe of points
df = pd.DataFrame([[0, 0, 0], [1, 2, 3], [2, 2, 2], [3, 2, 1] ], columns = list("XYZ"))
#your polygon points
polygon1_list = [(1, 1), (1, 3), (3, 3), (3, 1)]
#adding a column that contains a boolean variable for each point
df["polygon1"] = df.apply(lambda row: Polygon(polygon1_list).contains(Point(row["X"], row["Y"])), axis = 1)
print(df)
```
Ausgang für mein Spielzeug Datensatz
```
   X  Y  Z  polygon1
0  0  0  0   False
1  1  2  3   False
2  2  2  2    True
3  3  2  1   False
```
In formschönen, contains wirklich bedeutet, dass innerhalb des Polygons, das schließt die Grenze. Wenn Sie möchten, dass auch die Grenze, die Sie verwenden sollten intersects
```
df["polygon1"] = df.apply(lambda row: Polygon(polygon1_list).intersects(Point(row["X"], row["Y"])), axis = 1)
```
Nun die Antwort auf Ihre Frage ist einfach. Legen Sie einfach die Zeilen, die False in diese neue Spalte:
```
df = df.drop(df[~df["polygon1"]].index)
```
Leider haben Sie immer noch eine Schleife über die polygon-Liste. Es wäre interessant, wenn jemand wüsste einen Weg, wie um zu testen, alle Punkte und alle Polygone, die ohne eine (explizite) Schleife. Ich habe gesehen, eine MultiPolygon-Konstruktor der Klasse auf Ihrer website, also vielleicht eine Kombination aller Polygone in einer Klasse würde den trick tun. Aber testen Sie im Vorfeld, dass dies eine gültige Auswahl. Ein MultiPolygon ist ungültig, wenn seine Mitglieder berühren sich mit einer unendlichen Anzahl von Punkten entlang einer Linie.

Edit: Scheinbar ist in Python 2.7 funktioniert das nicht. Sehen akozi Antwort für einen 2.7 kompatibel zu beantworten.
- Das wirklich funktioniert der trick! Danke für die Anregung! Angesichts der point-of-Rutger Kassies ich bin jetzt auf der Suche nach Möglichkeiten, dies zu tun, vollständig mit geopandas. Vielleicht eine sauberere Lösung mit weniger Abhängigkeiten
- Ich habe gesehen, gestern und dachte an diesen thread. Vielleicht Hilfe für Sie. Viel Glück mit Ihrem Projekt.
- Vielen Dank nochmal, check it out
- Für jeden, der versucht dies in Python 2.7 müssen Sie contains_points statt contains. Es scheint auch Probleme mit einzelnen Elementen. Habe ich ein funktionsfähiges Beispiel bearbeitet diese Antwort in eine Antwort auf meine eigenen.
InformationsquelleAutor Mr. T
1

Ich hatte Mühe, imitiert die genau die Lösung, Mr. T vorgeschlagen, in Python 2.7. So, hier ist der kleine Unterschied, den ich machen musste, um es arbeiten in Python 2.7.
```
from shaply.geometry.polygon import Polygon
inside = Polygon(poly_points).contains_points(zip(df.X.values, df.Y.values))
df['inside'] = inside
df = df.drop(df[~df['inside']].index)
```
Scheint es, dass die alte version von contains_points hatte Probleme bei der Ausführung mit einem einzigen Punkt. So habe ich es bis zu Lesen Sie alle Punkte und hängen Sie diese Liste in einer neuen Spalte.

InformationsquelleAutor akozi

Schreibe einen Kommentar Antworten abbrechen