Mam zestaw ~ 36 000 wielokątów reprezentujących partycję (~ hrabstwa) kraju. Mój skrypt Pythona otrzymuje wiele punktów: pointId, longitude, latitude.

Dla każdego punktu chcę wysłać backInd, polygonId. Dla każdego punktu zapętlenie wszystkich wielokątów i użycie myPoint.within (myPolygon) jest dość nieefektywne.

Przypuszczam, że kształtna biblioteka oferuje lepszy sposób przygotowania wielokąta, aby znalezienie wielokąta dla punktu stało się ścieżką drzewa (kraj, region, podregion, ...)

Oto mój kod do tej pory:

import sys
import os
import json
import time
import string
import uuid

py_id = str(uuid.uuid4())

sys.stderr.write(py_id + '\n')
sys.stderr.write('point_in_polygon.py V131130a.\n')
sys.stderr.flush()

from shapely.geometry import Point
from shapely.geometry import Polygon
import sys
import json
import string
import uuid
import time

jsonpath='.\cantons.json'
jsonfile = json.loads(open(jsonpath).read())

def find(id, obj):
    results = []

    def _find_values(id, obj):
        try:
            for key, value in obj.iteritems():
                if key == id:
                    results.append(value)
                elif not isinstance(value, basestring):
                    _find_values(id, value)
        except AttributeError:
            pass

        try:
            for item in obj:
                if not isinstance(item, basestring):
                    _find_values(id, item)
        except TypeError:
            pass

    if not isinstance(obj, basestring):
        _find_values(id, obj)
    return results


#1-Load polygons from json  
r=find('rings',jsonfile)
len_r = len(r)

#2-Load attributes from json
a=find('attributes',jsonfile)

def insidePolygon(point,json):
        x=0
    while x < len_r :
            y=0
            while y < len(r[x]) :
        p=Polygon(r[x][y])
        if(point.within(p)):
            return a[y]['OBJECTID'],a[y]['NAME_5']
        y=y+1
        x=x+1
    return None,None


while True:
    line = sys.stdin.readline()
    if not line:
        break
    try:
        args, tobedropped = string.split(line, "\n", 2)

        #input: vehicleId, longitude, latitude
        vehicleId, longitude, latitude = string.split(args, "\t")

        point = Point(float(longitude), float(latitude))
        cantonId,cantonName = insidePolygon(point,r)

        #output: vehicleId, cantonId, cantonName
        # vehicleId will be 0 if not found
        # vehicleId will be < 0 in case of an exception
        if (cantonId == None):
            print "\t".join(["0", "", ""])
        else:
            print "\t".join([str(vehicleId), str(cantonId), str(cantonName)])

    except ValueError:
        print "\t".join(["-1", "", ""])
        sys.stderr.write(py_id + '\n')
        sys.stderr.write('ValueError in Python script\n')
        sys.stderr.write(line)
        sys.stderr.flush()
    except:
        sys.stderr.write(py_id + '\n')
        sys.stderr.write('Exception in Python script\n')
        sys.stderr.write(str(sys.exc_info()[0]) + '\n')
        sys.stderr.flush()
        print "\t".join(["-2", "", ""])