Suite

Utiliser la fonction dans la requête SQL avec ArcPy ?

Utiliser la fonction dans la requête SQL avec ArcPy ?


J'essaie de joindre des points à des polygones en fonction de la correspondance de leurs adresses et de la distance qui les sépare est inférieure à 2500 pieds. J'ai défini une fonction haversine qui calcule la distance entre le point et le polygone en utilisant leurs coordonnées. Cependant, j'obtiens une erreur indiquant que le SQL n'est pas valide.

import arcpy import harversine arcpy.env.overwriteOutput=True arcpy.env.workspace=r"C:mygdb.gdb" where="pointshp.CAdd = polyshp.P_Address et harversine.myhaversine(pointshp.NewLong,pointshp.NewLat,polyshp .polyX,polyshp.polyY) < 2500" keyfield="pointshp.OBJECTID" arcpy.MakeQueryTable_management(["pointshp","polyshp"],"joined","USE_KEY_FIELDS",keyfield,"",where) arcpy.CopyFeatures_management( "rejoint","rejoint")

L'erreur que j'obtiens :

Une instruction SQL non valide a été utilisée. [rejoint]

lehaversinefonction que j'utilise qui est enregistrée dans un script séparé appeléharversine.pyc'est-à-dire que ce n'est pas une faute d'orthographe.

import math def myhaversine(lon1, lat1, lon2, lat2) : # convertir les degrés décimaux en radians lon1, lat1, lon2, lat2 = map(math.radians, [lon1, lat1, lon2, lat2]) # formule haversine dlon = lon2 - lon1 dlat = lat2 - lat1 a = math.sin(dlat/2)**2 + math.cos(lat1) * math.cos(lat2) * math.sin(dlon/2)**2 c = 2 * math.asin(math.sqrt(a)) r = 20887680 # Rayon de la terre en kilomètres. Utilisez 3956 pour les miles, 20 887 680 pieds renvoient math.ceil (c * r)

Je pense que votre meilleur pari sera de déterminer d'abord toutes les adresses qui correspondent aux critères de distance, puis de filtrer celles-ci. Voici une façon de procéder. Déterminez les valeurs des points X et Y à partir des champs et ajoutez-les à un dictionnaire avec la clé étant l'adresse. Ensuite, effectuez les calculs par adresse tout en itérant à travers les polygones, et si l'adresse répond aux critères de distance, ajoutez l'adresse à une liste de « bonnes adresses ». Enfin, incorporez ces adresses dans la clause where. Tout est tabulaire, donc cela devrait être relativement rapide.

Le code (non testé) :

import arcpy import harversine arcpy.env.overwriteOutput = True arcpy.env.workspace = r"C:mygdb.gdb" print "création du dictionnaire de points" #Dictionnaire avec XYs par adresse à partir du point pntDi = dict ([(addr, (x , y)) pour addr, x, y dans arcpy.da.SearchCursor ("pointshp", ["CAd", "NewLong", "NewLat"]) si addr et x et y]) #liste vide pour les bonnes adresses goodAddresses = [] #iterate polygon table print "itération des polygones et calcul" avec arcpy.da.SearchCursor ("polyshp", ["P_Address", "polyX", "polyY"]) comme curseur : pour addr, polyX, polyY dans le curseur : si pas addr ou pas polyX ou pas polyY : continuer si pas addr dans pntDi : continuer pntX = pntDi[addr][0] pntY = pntDi[addr][1] #Calculer la valeur value = harversine.myhaversine(pntX, pntY, polyX, polyY) #Ajouter une bonne valeur à la liste d'adresses si valeur < 2500 : goodAddresses += [addr] #Nouvelle clause whereStr = "', '".join (goodAddresses) where = "pointshp.CAd = polyshp.P_Address et pointshp. CAjouter dans ('{0}')".format (addressStr)