Point dans le polygone avec geoJSON en Python
J'ai une base de données geoJSON avec beaucoup de polygones (secteurs de recensement spécifiquement) et j'ai beaucoup de points longs et lat.
J'espère qu'il existerait un code Python efficace pour identifier le secteur de recensement dans lequel se trouve une coordonnée donnée, mais jusqu'à présent, mon googling n'a rien révélé.
Merci!
3 réponses
J'ai trouvé un article intéressant décrivant comment faire exactement ce que vous cherchez à faire.
TL;DR: Utiliser Galbées
, Vous trouverez ce code à la fin de l'article:
import json
from shapely.geometry import shape, Point
# depending on your version, use: from shapely.geometry import shape, Point
# load GeoJSON file containing sectors
with open('sectors.json') as f:
js = json.load(f)
# construct point based on lon/lat returned by geocoder
point = Point(-122.7924463, 45.4519896)
# check each polygon to see if it contains the point
for feature in js['features']:
polygon = shape(feature['geometry'])
if polygon.contains(point):
print 'Found containing polygon:', feature
Une excellente option pour travailler avec ces types de données est PostGIS , Un extendeur de base de données spatiale pour PostgreSQL. Personnellement, je garde toutes mes données geo dans une base de données PostGIS, puis je les référence en python en utilisant psycopg2. Je sais que ce n'est pas Python pur, mais il a des avantages de performance incroyables (discutés ci-dessous) sur Python pur.
PostGIS a une fonctionnalité intégrée pour déterminer si un point ou une forme se trouve dans une autre forme. La bonne documentation sur le la fonction ST_Within développe cet exemple simple:
SELECT
ST_WITHIN({YOUR_POINT},boundary)
FROM census;
-- returns true or false for each of your tracts
Le bénéfice que vous gagnerez de PostGIS que vous n'obtiendrez probablement pas ailleurs est l'indexation, ce qui peut améliorer votre vitesse 1 000 x [1], ce qui le rend meilleur que même le meilleur programme c écrit (à moins que le programme C crée également un index pour vos données). La base de données, lorsqu'elle est correctement configurée, met en cache des informations sur vos tracts, et lorsque vous demandez si un point se trouve dans un tract, il n'aura pas à tout rechercher... il peut profitez de son index.
Obtenir des données dans et hors de PostGRES est assez simple. Un excellent tutoriel qui vous guidera à travers les bases de PostGIS avec des exemples de jeux de données pas trop différents du vôtre peut être trouvé ici. C'est assez long, mais si vous êtes nouveau à PostGIS (comme je l'étais), vous serez très diverti et excité tout le temps:
Http://workshops.boundlessgeo.com/postgis-intro/
[1] L'indexation a diminué une recherche de voisin le plus proche dans l'un de mes bases de données énormes (20 m de 53 secondes à 8,2 millisecondes.
On ne peut pas avoir de code géométrique très rapide en Python. Au lieu de cela, l'approche habituelle consiste à utiliser une bibliothèque C/C++ Rapide avec des wrappers Python.
Par exemple, vous pouvez commencer avec CGAL - une bibliothèque géométrique c++ très complète. Il a des liaisons Python pour la plupart de ses routines, voir le lien http://code.google.com/p/cgal-bindings/.