Many Points to Few...

Discussion created by Cjkinzett on Dec 9, 2013
Latest reply on Dec 11, 2013 by Cjkinzett
Hello, i am unsure if this is the correct place for this but it is an arcpy problem and i am new here, feel free to redirect me.

Anyway I am trying to reduce a large dataset of about 4.6 Million points down to a few hundred thousand by removing duplicates. I have been running this in arcpy as its more stable and faster, and who doesn't love some programming.

So far my script is less than optimal but essentially it;

Clips points by New Zealand regions
     clips points by mesh-blocks within this region
     deletes identical points within these mesh-blocks by buffering 20 metres out from the '*shape'
     appends it to a list of new points.

It works quite well and am about halfway however it aborts at the same point everytime,

11.95% Complete

i have tried repair geometry for all datasets and have tried changing from meshblocks to area-units but to no avail.Here is most of the code. Cheers in advance!

"""Run check geometry on regions"""
gp.addmessage("Checking Geometries")
arcpy.CheckGeometry_management ([NZAreaUnits,TAs,"points"], "ManyPointsTable")
count = arcpy.GetCount_management("ManyPointsTable")
gp.addmessage(str(count) + " errors found")
if count >= 1:
    gp.addmessage("...Repairing ")
    arcpy.RepairGeometry_management (NZAreaUnits)
    arcpy.RepairGeometry_management (TAs)
    arcpy.RepairGeometry_management ("Points")
    gp.addmessage("Search Cursor")
"""Create Search cursor for Large Regions"""
seaCur = arcpy.SearchCursor(TAs)
region = 8
number = float(arcpy.GetCount_management(TAs).getOutput(0))

"""Make Small Features Layer and Select if they contain points then clip and use for Small regions"""
gp.addmessage("Making Feature Layer with only valued Features")
arcpy.MakeFeatureLayer_management(NZAreaUnits, "areunits")
arcpy.SelectLayerByLocation_management ( "areunits", "WITHIN_A_DISTANCE", "points" , "20 METERS" )
arcpy.FeatureClassToFeatureClass_conversion ("areunits", "C:\\Users\\Cody\\Documents\\GISData\\AvlSH2013Take2.gdb", "ClippedAU")
NZAreaUnits ="ClippedAU"

#Regional Clip
for Regions in seaCur:
     percentcomplete = float((region/number)* 100)
     region += 1
     #Clip the features by row
     geom = Regions.shape
     gp.addmessage("Running Clip on -" + Regions.getValue("NAME"))
     arcpy.Clip_analysis("points" , geom , "Clip")
     arcpy.Clip_analysis(NZAreaUnits , geom , "AUClip")
     gp.addmessage("Starting Region" + " "+str(region))
     SearchCur = arcpy.SearchCursor("AUClip")
     #Create New Dataset and name
     outset = "Finalset" + str(region)
     arcpy.CreateFeatureclass_management("C:\\Users\\Cody\\Documents\\GISData\\AvlSH2013Take2.gdb", outset, "POINT", "points", "DISABLED", "DISABLED", "points")
     number1 = float(arcpy.GetCount_management("AUClip").getOutput(0))
     gp.addmessage(str("{0:.2f}".format(round(percentcomplete,2))) + "% Total Complete")
     AU = 0
     #AU Clip
     for AreUnit in SearchCur:
         AU += 1
         percomplete = float((AU/number1)* 100)       
         geom = AreUnit.shape
         arcpy.Clip_analysis("Clip" , geom , "Sclip") 
         #Delete Identical
         arcpy.MakeFeatureLayer_management("Sclip", "Layer")
         arcpy.DeleteIdentical_management("Layer", "Shape", "20 METERS")
         #Merge Data
         arcpy.Merge_management(["Layer", outset ], "Merge")
         arcpy.CopyFeatures_management ("Merge", outset)
         gp.addmessage(str("{0:.2f}".format(round(percomplete,2))) + "% Complete")