Community

Need help to facilitate the process

384
2
03-26-2012 03:47 PM
zhengniu1
by
New Contributor
‎03-26-2012 03:47 PM
Hi, python Supermen, I got a pseudo code as described below. I am trying to make it works, but running is really slow. Please help to facilitate it if you dont mind.

The input table is as shown in attachment. 

C = 0 (set of rep. points = empty)
for sj in S: (All GPS points involved in all trajectories)
 find points are within (4m) to sj & are  not represented by any rep. in C
 calculate the centroid c*1 of these points
 if no points are within (4m) to sj, then sj will represent itself (or marked as outlier/noise to be removed later)
 find points (k) are within (2m) to c*1
 for sk in k
  if sk is not represented
  assign sk to c*1
  if sk is represented by c*2 
  compare distances  dk1 & dk2
       if dk1 < dk2 
   assign sk to c*1 
        If all neighbors of sk are represented, then randomly choose sj from the remaining un-represented points
 
               Iterate until all sj are represented 

C stores all rep. point (structure: C:[Id, Pc_X, Pc_Y, dir, ncp:{all pts within 2m (less or equal to) of  c*i}] s) *
Tags (2)
Preview file
222 KB
0 Kudos
2 Replies
zhengniu1
by
New Contributor
‎03-26-2012 03:51 PM
Part 1

import arcpy
import math
import time 

arcpy.env.workspace=r'C:\TempArcGISCalculate\Scratch.gdb'
arcpy.env.overwriteOutput = True
t0=time.clock()
OriginNearTable="Section01Smoothed2_GenerateN"
NearTable = "TableView_Section01Smoothed2_GenerateN"
arcpy.MakeTableView_management(OriginNearTable, NearTable)
result = "Section01ClusterSample11"

arcpy.CreateFeatureclass_management("C:/TempArcGISCalculate/Scratch.gdb", result,"POINT","","","", "C:\TempArcGISCalculate\WGS_1984_UTM_Zone_17N.prj")
arcpy.AddField_management(result, "ClusterID","LONG")
arcpy.AddField_management(result, "C_LON_X","DOUBLE")
arcpy.AddField_management(result, "C_LAT_Y","DOUBLE")
arcpy.AddField_management(result, "C_Bearing","DOUBLE")
arcpy.AddField_management(result, "NearbyPtID","LONG")
arcpy.AddField_management(result, "NearbyLON_X","DOUBLE")
arcpy.AddField_management(result, "NearbyLAT_Y","DOUBLE")
arcpy.AddField_management(result, "NearbyBearing","DOUBLE")
ptCursor = arcpy.InsertCursor(result)

InFIDlst=[]
NearFIDlst=[]
NearXlst=[]
NearYlst=[]
NearBearinglst=[]
TargetXlst=[]
TargetYlst=[]
TargetBearinglst=[]
TargetRouteIDlst=[]
TargetTracelst=[]
Reduncelst =[]

Currcur = arcpy.SearchCursor(NearTable) # Current Cursor for field: IN_FID
Prevcur = arcpy.SearchCursor(NearTable) # Previous Cursor for field: IN_FID
ccur = Currcur.next()
pcur = Prevcur.next()

prevINFID = ccur.getValue('IN_FID')
InFIDlst.append(ccur.getValue('IN_FID'))
NearFIDlst.append(ccur.getValue('NEAR_FID'))
NearXlst.append(ccur.getValue('NEAR_X'))
NearYlst.append(ccur.getValue('NEAR_Y'))
NearBearinglst.append(ccur.getValue('Bearing_1'))
TargetXlst.append(ccur.getValue('LON_X'))
TargetYlst.append(ccur.getValue('LAT_Y'))
TargetBearinglst.append(ccur.getValue('Bearing'))
TargetRouteIDlst.append(ccur.getValue('RouteId'))
TargetTracelst.append(ccur.getValue('Trace'))

ccur=Currcur.next()

flag=0
IDcount = 1 # Cluster ID
Rcount = 0 # To count number of target point which is already to be neighoubr of other points
Scount = 0 #

while pcur:
    if ccur:
        currentINFID = ccur.getValue('IN_FID')
    else:#When it reached the end of the file
        currentINFID = -1
    if  prevINFID not in Reduncelst and (currentINFID!= prevINFID or not ccur):
        Clustlst=[]
        n = len(InFIDlst)                
        Clustlst.append((pcur.getValue('IN_FID'),pcur.getValue('LON_X'),pcur.getValue('LAT_Y'),pcur.getValue('Bearing'),pcur.getValue('IN_FID'),pcur.getValue('LON_X'),pcur.getValue('LAT_Y'),pcur.getValue('Bearing')))
        if n<>0: # avoid Clustlst = [] meaning the neighbor points of the target point in IN_FID are used by other target points
            for i in range(n):
                Clustlst.append((InFIDlst,TargetXlst,TargetYlst,TargetBearinglst,NearFIDlst, NearXlst,NearYlst,NearBearinglst))
                if NearFIDlst not in Reduncelst:
                    Reduncelst.append(NearFIDlst)
            TX = 0.0
            TY = 0.0
            TDir = 0.0
            CX = 0.0
            CY = 0.0
            CDir = 0.0
            for j in range(len(Clustlst)):                
                TX = TX + Clustlst[5]
                TY = TY + Clustlst[6]
                TDir = TDir + Clustlst[7]
            CX=TX/len(Clustlst)
            CY=TY/len(Clustlst)
            CDir =TDir/len(Clustlst)
     #Make a in-memory table to store all the points related to the point in current cluster
            arcpy.SelectLayerByAttribute_management(NearTable, "NEW_SELECTION",  "IN_FID in (" + str(NearFIDlst)[1:-1] + ")") # when NearFIDlst !=[]
            arcpy.CopyRows_management(NearTable, "in_memory/selectedpoints")
            tempcur = arcpy.SearchCursor("in_memory/selectedpoints")
    #See if there is any points not included in current cluster
            for tcur in tempcur: # Search neighbor point's neighbor
                if tcur.NEAR_FID in NearFIDlst or tcur.NEAR_FID==InFIDlst[0]: # if the neighbor of the target's neighbor point is already included or the neighbor of the target's neighor point is the target itself
                    continue    # ignore it and then continue
                elif ((tcur.NEAR_X-CX)**2+(tcur.NEAR_Y-CY)**2)<= 16: # if the neighbor of the target's neighbor point is within the 4-m buffer of the preliminary cluster center
                    comparecur = arcpy.SearchCursor(result)
                    for compcur in comparecur: # Search classified points represented by the generated cluster center
                        #When {pi} is found in the "Section01ClusterSample"
                        if compcur.NearbyPtID == tcur.NEAR_FID: # if the neighbor of the target's neighbor point is already represented by one final cluster center
                            # Compare distances from the neighbor of the target's neighbor point to that final cluster center and the preliminary cluster center
                            dist2 = math.sqrt((compcur.getValue('NearbyLON_X') - compcur.getValue('C_LON_X'))**2+(compcur.getValue('NearbyLAT_Y') - compcur.getValue('C_LAT_Y'))**2)
                            dist1 = math.sqrt((compcur.getValue('NearbyLON_X')-CX)**2+(compcur.getValue('NearbyLAT_Y')-CY)**2)    
                            if dist1<dist2:
                                arcpy.MakeTableView_management(result,"tempViewSection01ClusSample")
                                arcpy.SelectLayerByAttribute_management("tempViewSection01ClusSample","NEW_SELECTION",'"NearbyPtID" =%i' %compcur.NearbyPtID)
                                arcpy.DeleteRows_management("tempViewSection01ClusSample")
                                flag = 1 #flag state 1 means {pi} closer to current siling (--- preliminary cluster center)
                                break
                            else:
                                flag = 2 #flag state 2 means {pi} closer to original siling (--- the final cluster center)
                                break
                        else:
                            flag = 0 # if no compcur.NearbyPtID qualified as tcur.NEAR_FID (zhao ren, mei zhao dao, meizhaodao, meizhaodao...), go to step: if flag<>2
                            #flag state 0 means {pi} not presented by any siling
                    if flag<>2:
                        Clustlst.append((Clustlst[0][0], Clustlst[0][1],Clustlst[0][2],Clustlst[0][3],tcur.NEAR_FID,tcur.NEAR_X,tcur.NEAR_Y,tcur.Bearing_1))
                        if tcur.NEAR_FID not in Reduncelst:
                            Reduncelst.append(tcur.NEAR_FID)
0 Kudos
zhengniu1
by
New Contributor
‎03-26-2012 03:51 PM
Part 2

        else:
            CX = Clustlst[0][5]
            CY = Clustlst[0][6]
            CDir = Clustlst[0][7]
            print "Scount:", Scount, "ClusterID",IDcount, "X", CX, "Y",CY
            Scount=Scount+1
        k = 0
        while k <len(Clustlst):
            row=ptCursor.newRow()
            point = arcpy.CreateObject('POINT')
            point.X=CX
            point.Y=CY
            row.Shape=point
            row.ClusterID = IDcount
            row.C_LON_X = CX
            row.C_LAT_Y = CY
            row.C_Bearing = CDir
            row.NearbyPtID = Clustlst[4]
            row.NearbyLON_X = Clustlst[5]
            row.NearbyLAT_Y = Clustlst[6]
            row.NearbyBearing = Clustlst[7]
            ptCursor.insertRow(row)
            k=k+1          
        InFIDlst=[]
        NearFIDlst=[]
        NearXlst=[]
        NearYlst=[]
        NearBearinglst=[]
        TargetXlst=[]
        TargetYlst=[]
        TargetBearinglst=[]
        TargetRouteIDlst=[]
        TargetTracelst=[]
        IDcount = IDcount+1
    elif currentINFID!= prevINFID or not ccur:
        InFIDlst=[]
        NearFIDlst=[]
        NearXlst=[]
        NearYlst=[]
        NearBearinglst=[]
        TargetXlst=[]
        TargetYlst=[]
        TargetBearinglst=[]
        TargetRouteIDlst=[]
        TargetTracelst=[]
        Rcount = Rcount+1
    if ccur:
        tmp1=ccur.getValue('NEAR_FID')
        tmp2=ccur.getValue('IN_FID')
        if tmp1 in Reduncelst and tmp2 not in Reduncelst:#To compare before adding a point and avoid centers in "Rcount" Siling
            comparecur = arcpy.SearchCursor(result)
            for compcur in comparecur:
                if compcur.NearbyPtID == ccur.NEAR_FID:
                    dist2 = math.sqrt((compcur.getValue('NearbyLON_X') - compcur.getValue('C_LON_X'))**2+(compcur.getValue('NearbyLAT_Y') - compcur.getValue('C_LAT_Y'))**2)
                    dist1 = math.sqrt((compcur.getValue('NearbyLON_X')-ccur.getValue('LON_X'))**2+(compcur.getValue('NearbyLAT_Y')-ccur.getValue('LAT_Y'))**2)      
                    if dist1<dist2:
                        arcpy.MakeTableView_management(result,"tempViewSection01ClusSample")
                        arcpy.SelectLayerByAttribute_management("tempViewSection01ClusSample","NEW_SELECTION",'"NearbyPtID" =%i' %compcur.NearbyPtID)
                        arcpy.DeleteRows_management("tempViewSection01ClusSample")
                        InFIDlst.append(ccur.getValue('IN_FID'))
                        NearFIDlst.append(ccur.getValue('NEAR_FID'))
                        NearXlst.append(ccur.getValue('NEAR_X'))
                        NearYlst.append(ccur.getValue('NEAR_Y'))
                        NearBearinglst.append(ccur.getValue('Bearing_1'))
                        TargetXlst.append(ccur.getValue('LON_X'))
                        TargetYlst.append(ccur.getValue('LAT_Y'))
                        TargetBearinglst.append(ccur.getValue('Bearing'))
                        TargetRouteIDlst.append(ccur.getValue('RouteId'))
                        TargetTracelst.append(ccur.getValue('Trace'))
        elif tmp2 not in Reduncelst: #avoid adding if this is a center in "Rcount" to improve performance.
            InFIDlst.append(ccur.getValue('IN_FID'))
            NearFIDlst.append(ccur.getValue('NEAR_FID'))
            NearXlst.append(ccur.getValue('NEAR_X'))
            NearYlst.append(ccur.getValue('NEAR_Y'))
            NearBearinglst.append(ccur.getValue('Bearing_1'))
            TargetXlst.append(ccur.getValue('LON_X'))
            TargetYlst.append(ccur.getValue('LAT_Y'))
            TargetBearinglst.append(ccur.getValue('Bearing'))
            TargetRouteIDlst.append(ccur.getValue('RouteId'))
            TargetTracelst.append(ccur.getValue('Trace'))
    prevINFID = currentINFID 
    ccur = Currcur.next()
    pcur = Prevcur.next()
del ccur,Currcur, pcur,Prevcur,row, ptCursor,tcur, tempcur,compcur, comparecur
print "Elapsed Time:", time.clock() - t0, "Seconds", "IDcount=", IDcount, "Rcount=", Rcount, "Scount=", Scount
0 Kudos