GeoPandas intersects ne trouve aucune intersection

Sep 30 2020

Je veux vérifier si la couche de points que j'ai intersecte avec la couche de polygones que j'ai, en tant que colonne booléenne dans le dataframe de points.
J'ai deux dataframes GeoPandas, le premier est de nombreux points et ressemble à ceci:

>>>ID  geometry
0  12  POINT (5.0279 7.4547)
1  45  POINT (6.6539 12.139)
...

et le deuxième dataframe est une couche de nombreux polygones différents qui ressemble à ceci:

>>>name     code   geometry
0  Desert   12     POLYGON ((5.52013 13.8902, 5.5265 13.892,...)
1  Water    24     POLYGON ((5.53756 13.88472, 5.5291 13.8791,...)
...

J'essaie de vérifier s'il y a une intersection entre la couche de points et la couche de régions. Pour cela, j'ai déterminé les crs puis utilisé des intersections:

regions=regions.to_crs({'init': 'epsg:4326'})
points=points.set_crs({'init': 'epsg:4326'})
inter=points.geometry.intersects(regions.geometry)

A Le script s'exécute avec les avertissements suivants:

FutureWarning: '+ init =:' syntax is deprecated. ':' is the preferred initialization method. When making the change, be mindful of axis order changes: https://pyproj4.github.io/pyproj/stable/gotchas.html#axis-order-changes-in-proj-6 return _prepare_from_string(" ".join(pjargs)) /opt/conda/lib/python3.8/site-packages/geopandas/base.py:39: UserWarning: The indices of the two GeoSeries are different.
warn("The indices of the two GeoSeries are different.")

Then when I check the results the only value is False, like all th epoints do not intersect:

inter.unique().tolist()
>>>[False]

* J'ai vu sur QGIS qu'il y a des points qui se croisent et qu'il y a des points qui ne le font pas donc il n'y a aucun moyen que ce résultat soit vrai

* J'ai vérifié les dtypes - chacun de mes geodataframes a une colonne qui est la géométrie et appelée géométrie.

Mon objectif final: ajouter une nouvelle colonne dans le geodataframe des points qui dira si elle intersecte les régions ou non.

Réponses

4 martinfleis Sep 30 2020 at 15:13

Intersects est une opération par ligne qui aligne les deux GeoSeries et check est une géométrie A dans la ligne 0 coupe la géométrie B dans la ligne 0.

Pour vérifier s'il y a une intersection, utilisez l'index spatial.

import numpy as np

inp, res = regions.sindex.query_bulk(points.geometry, predicate='intersects')
points['intersects'] = np.isin(np.arange(0, len(points)), inp)

Concernant l'avertissement CRS, supprimez simplement une initpartie et passez-la sous forme de chaîne.

regions=regions.to_crs('epsg:4326')
points=points.set_crs('epsg:4326')

query_bulknécessite geopandas 0.8 et n'est pas encore documenté, mais encapsule pygeos.STRtree.query_bulk, vous pouvez donc vérifier que la documentation -https://pygeos.readthedocs.io/en/latest/strtree.html#pygeos.strtree.STRtree.query_bulk.