Python Near Analysis (No Advanced Licence)
I need to determine the shortest distance between each point within a single feature class. Near Analysis would have worked perfectly, but unfortunately you need to have an Advanced ArcGIS Licence and I only have a Standard Licence.
I found the following code on Stack Exchange written by PolyGeo:
import arcpy,math
# Set variables for input point feature classes and output table
ptFC1 = "C:/temp/test.gdb/PointFC1"
ptFC2 = "C:/temp/test.gdb/PointFC2"
outGDB = "C:/temp/test.gdb"
outTableName = "outTable"
outTable = outGDB + "/" + outTableName
arcpy.env.overwriteOutput = True
# Create empty output table
arcpy.CreateTable_management(outGDB,outTableName)
arcpy.AddField_management(outTable,"INPUT_FID","LONG")
arcpy.AddField_management(outTable,"NEAR_FID","LONG")
arcpy.AddField_management(outTable,"DISTANCE","DOUBLE")
# Create and populate two dictionaries with X and Y coordinates for each
# OBJECTID in second feature class using a SearchCursor
ptFC2XCoordDict = {}
ptFC2YCoordDict = {}
with arcpy.da.SearchCursor(ptFC2,["OBJECTID","SHAPE@XY"]) as cursor:
for row in cursor:
ptFC2XCoordDict[row[0]] = row[1][0]
ptFC2YCoordDict[row[0]] = row[1][1]
# Open an InsertCursor ready to have rows written for each pair of OBJECTIDs
iCursor = arcpy.da.InsertCursor(outTable,["INPUT_FID","NEAR_FID","DISTANCE"])
# Use a SearchCursor to read the rows (and X,Y coordinates) of the first
# feature class
with arcpy.da.SearchCursor(ptFC1,["OBJECTID","SHAPE@XY"]) as cursor:
for row in cursor:
x1 = row[1][0]
y1 = row[1][1]
for i in range(len(ptFC2XCoordDict)):
x2 = ptFC2XCoordDict[i+1]
y2 = ptFC2YCoordDict[i+1]
# Prepare and insert the InsertCursor row
iRow = [row[0],i+1,math.sqrt((x2-x1)*(x2-x1) + (y2-y1)*(y2-y1))]
iCursor.insertRow(iRow)
del iCursor
print "Done!"The following reports the distance of the first point to the rest of the points:
I would like to amend the following to iterate through each point within the feature class to report the shortest distance between the points, greater than zero of course. Any advice in amending the following will be appreciated. I'm trying to use the results from this as a way of validating that none of the points are within 3m or less from each other and if there are the filtered list is processed further to move the points away from each other within bounding box for each point.
Peter
distance from the first point in a point file to all other points in the file.
the rest of the problem is homework :
"""
Script: distance_einsum.py
Author:
Dan.Patterson@carleton.ca
Purpose: Yes
References:
Numpy Repository stuff
"""
import numpy as np
import arcpy
np.set_printoptions(edgeitems=3,linewidth=80,precision=2,suppress=True,threshold=10)
def einsum_dist(origin,destin):
"""using einsum"""
deltas = destin - origin
return np.sqrt(np.einsum('ij,ij->i', deltas, deltas))
def demo(a):
""" demo of one origin to all other points"""
orig = a['Shape'][0]
dest = a['Shape'][1:]
dist = einsum_dist(orig,dest)
out = np.zeros((len(dist)+1,3))
out[:,0] = a['Shape'][:,0] # the Xs
out[:,1] = a['Shape'][:,1] # the Ys
out[1:,2] = dist
return out
if __name__=="__main__":
"""sample run for distance demo"""
f = r"F:\A2_in_progress\NearDemo\Shapefiles\RandomPnts.shp"
shp_field = "Shape"
a = arcpy.da.FeatureClassToNumPyArray(f, shp_field)
b = demo(a)
# send it to a featureclass if desired
print("X, Y, distance from the first point \n{}".format(b))
result
>>> X, Y, distance from the first point
[[ 342059. 5022985. 0. ]
[ 342431. 5025573. 2614.6 ]
[ 341341. 5024714. 1872.16]
...,
[ 343020. 5022952. 961.57]
[ 343198. 5025388. 2659.27]
[ 342635. 5023786. 986.6 ]]
Now where are those close points and how many of them are there? slice and dice. b[:,2] is the distance column ... b[:,2] < 1000 are all the distances < 1000 m ... Therefore:
b[b[:,2]<1000] must be all the records where the distance is less than 1000m .. now what? save it to a featureclass. And yes you can get all fancy with inputs, retaining variable, etc. but that should suffice. It is pretty quick up about 50,000,000 origin destinations.
>>> b[b[:,2]<1000]
array([[ 342059. , 5022985. , 0. ],
[ 341427. , 5023380. , 745.28],
[ 342494. , 5022158. , 934.43],
...,
[ 342195. , 5023291. , 334.86],
[ 343020. , 5022952. , 961.57],
[ 342635. , 5023786. , 986.6 ]])
>>> len(b[b[:,2]<1000])
24
>>>
Check Numpy Repository or contact me for further information.
Thanks Dan, much appreciated. Will get back to you once I've done my homework 
HI Dan
Thanks, will give me a chance to figure it out first before you post the answer . Thanks for being such a great teacher. I have learnt a great deal from you already. My Python skill set has grown tremendously over the last 3 months.
