ArcPy find points that are paired in time and space

Question

I have a shapefile that that is time stamped and tracks two animals. I would like to create a field that denotes "together" based on proximity and location. And example of the data is:

time       animal        data
13:00:04      A            x
13:00:02      B            x
22:07:03      A            x
23:08:04      B            x

Together is defined as being within 20 meters of each other and points falling within 5 seconds of each other. So the first pair would be classed as 'together' but the second pair would not, assuming all these points were in the same location. I have looked for similar questions that use ArcPy but could not find any.

JamesLeversha · Answer 1 · 2016-06-30T22:35:21.990

Using arcpy you should be able to do this pretty easily, the distance is just a plane distance though not sure if that matters. I have also never used a cursor with itertools.combinations, but it works with the test dictionary that I created.

An answer here suggests you can

https://arcpy.wordpress.com/2012/02/01/find-overlapping-features/

Even if you can't use combinations on a cursor its definitely worth looking into, more info here:

https://docs.python.org/2.7/library/itertools.html

import itertools
from datetime import datetime
import math
#import arcpy

#check distance on plane
def dist(A,B):
    ax, ay = A
    bx, by = B
    return math.hypot(bx-ax, by-ay)

#build record dictionary
records ={1: ('A', datetime.strptime("13:00:04","%H:%M:%S"),(0,0)),
          2: ('B', datetime.strptime("13:00:02","%H:%M:%S"),(0,1)),
          3: ('A', datetime.strptime("22:07:03","%H:%M:%S"),(50,50)),
          4: ('B', datetime.strptime("23:08:04","%H:%M:%S"),(20,20)),
          5: ('C', datetime.strptime("23:08:01","%H:%M:%S"),(115,115)),
          6: ('C', datetime.strptime("13:00:01","%H:%M:%S"),(5,5))}

#this blanked section to use arcpy to get data instead of dicitonary
#fc = 'test.shp'
#fields = ['ANIMAL','TIME','SHAPE@XY', 'OTHERDATA']

#with arcpy.da.SearchCursor(fc,fields) as cursor:
    #for value1,value2 in itertools.combinations(cursor,2):
    #if combinations doesn't work here, you can create your own dict
    #counter = 0
    #for row in cursor:
        #counter +=1
        #records[counter] = row


#blank this for loop if you use search cursor section above
# this checks one record to another without checking itself
# or checking the reverse way i.e. this will check A,B A,C, B,C
# not A,B A,C, B,A B,C C,B
for record1, record2 in itertools.combinations(records,2):
        value1 = records[record1]
        value2 = records[record2]
        #check to make sure you arent checking the same animal
        if value1[0] != value2[0]:
            #get time difference
            tdelta = (value1[1] - value2[1]).total_seconds()
            #check that time difference is less than 5 seconds
            if abs(tdelta) <= 5:
                #get plane distance difference
                distDif = dist(value1[2],value2[2])
                #check if less than 20m on plane
                if distDif <= 20:
                    print 'CLOSE ANIMALS:'
                    print value1,value2, distDif, tdelta

UPDATE: Combinations does work on a da search cursor

import itertools
from datetime import datetime
import math
import arcpy

#check distance on plane
def dist(A,B):
    ax, ay = A
    bx, by = B
    return math.hypot(bx-ax, by-ay)   

#this blanked section to use arcpy to get data instead of dicitonary
fc = 'test.shp'
fields = ['ANIMAL','TIME','SHAPE@XY', 'OTHERDATA']

with arcpy.da.SearchCursor(fc,fields) as cursor:
    for value1,value2 in itertools.combinations(cursor,2):            
        #check to make sure you arent checking the same animal
        if value1[0] != value2[0]:
            #get time difference
            tdelta = (value1[1] - value2[1]).total_seconds()
            #check that time difference is less than 5 seconds
            if abs(tdelta) <= 5:
                #get plane distance difference
                distDif = dist(value1[2],value2[2])
                #check if less than 20m on plane
                if distDif <= 20:
                    print 'CLOSE ANIMALS:'
                    print value1,value2, distDif, delta

UPDATE: This should work, however I was unsure how you want to handle multiple close animals to each other so I just created a list for them nested in the dictionary. This will leave you with a long string in the match field that looks like '[(animal,timeDif,distdif),(animal2,timeDif2,distdif2)]'. Look at the ast module to convert this back into a list if you run a search cursor on that. Otherwise look at created multiple fields for each separate close animal.

import itertools
from datetime import datetime
import math
import arcpy

#check distance on plane
def dist(A,B):
    ax, ay = A
    bx, by = B
    return math.hypot(bx-ax, by-ay)   

#this blanked section to use arcpy to get data instead of dicitonary
fc = 'test.shp'
fields = ['ANIMAL','TIME','SHAPE@XY', 'OID@','OTHERDATA']

closeDict = {}
with arcpy.da.SearchCursor(fc,fields) as cursor:
    for value1,value2 in itertools.combinations(cursor,2):            
        #check to make sure you arent checking the same animal
        if value1[0] != value2[0]:
            #get time difference
            tdelta = (value1[1] - value2[1]).total_seconds()
            #check that time difference is less than 5 seconds
            if abs(tdelta) <= 5:
                #get plane distance difference
                distDif = dist(value1[2],value2[2])
                #check if less than 20m on plane
                if distDif <= 20:

                    closeVal1 = closeDict.get(value1[3],[])
                    closeAnimal1 = value2[0],distDif,tdelta
                    closeVal1.append(closeAnimal1)
                    closeDict[value1[3]] = closeVal1

                    closeVal2 = closeDict.get(value2[3],[])
                    closeAnimal2 = value1[0],distDif,tdelta
                    closeVal2.append(closeAnimal2)
                    closeDict[value2[3]] = closeVal2

fields = ['OID@','MATCHFIELD']                    
with arcpy.da.UpdateCursor(fc,fields) as cursor:
    for row in cursor:
        try:
            row[1] = str(closeDict[row[0]])
            cursor.updateRow(row)

        except KeyError:
            pass

Thanks for this. I haven't been able to try it out yet, I've been caught up in another project. I will report back when I give it a shot! Thanks again! — Yolo_chicken, Jun 27 '16 at 12:56
How do I save this output rather than print it in the ArcPy module within ArcGIS? — Yolo_chicken, Jun 30 '16 at 12:18
@MB_analyst What do you want to save the results as, check out the answer by dslamb, it contains an update cursor which is useful. — JamesLeversha, Jun 30 '16 at 12:21
Ideally, results would add another field denoting "close" or not. To allow further analysis of when animals were within proximity — Yolo_chicken, Jun 30 '16 at 14:53

score 2 · Answer 2 · edited Apr 13 '17 at 12:33

2

Sort table by time, if it is not already a case and assign good unique IDs
Convert points to line using field animal as line ID.
Split lines at points
Remove lines longer than 20 m
Assign FROM node and TO node IDs to segment using something like this
Transfer time at FROM node and TO node to table of segments
Find difference in seconds. Expect trouble at midnight

edited Apr 13 '17 at 12:33

Community

1

answered Jun 21 '16 at 21:10

FelixIP

22,922
3
29
61

When I split lines at points, it's only giving me a few lines. Any thoughts on this? – Yolo_chicken Jul 11 '16 at 14:50
I'd do everything using shape files, but save time to text field before. For n points you should get n-1 lines. If this isn't a case, increase tolerance or whatever in split tool – FelixIP Jul 11 '16 at 19:46

dslamb · Answer 3 · 2016-06-23T12:45:05.317

You could follow this basic procedure. First, I would move your point feature class from a shapefile to a feature class in a file geodatabase. Then I would make your time column a datetime column so you can sort your time column in ascending order. This will make the process easier. Also, assign a unique ID to each point so you can match your pairs based on unique ID. If they are positive, then you can use -1 to flag a point that does not have a pair. Using the da module should return the time as datetime objects.

Assuming you have a unique ID field, time field, and match field.

import arcpy

tempList = []
resultDict = {}
with arcpy.da.SearchCursor(fc,["SHAPE@","UID","TIMEFIELD"]) as sc:
    for row in sc:
        tempList.append([row[1],row[0],row[2]]) #list of lists, id, pointgeometry, time

Then go through tempList and store the match pair in the dictionary

for indx, val in enumerate(tempList):
    cp = val[1]
    ct = val[2]
    matchID = -1
    if not val[0] in resultDict.keys():
        if indx < len(tempList):
            for item in tempList[indx:]:
                if item[0] != val[0]:
                    d = cp.distanceTo(item[1])
                    td = item[2] - ct  #should be a time delta object
                    if d <=20 and td.total_seconds() <=5:
                        matchID = item[0]
                        break
                    if td.total_seconds() > 5:
                        break #no need to continue looping since it exceeds five seconds
         resultDict[matchID]=val[0]
         resultDict[val[0]]=matchID #assuming that the matches are equal for each point

Then write the results:

with arcpy.da.UpdateCursor(fc,["UID","MATCHFIELD"]) as uc:
    for row in uc:
        row[1] = resultDict[row[0]]
        uc.updateRow(row)

I realize this may not be the most efficient path to take. I'm not familiar with the itertools approach, and that may be better than the enumerate function.

Edited the code to add a break if the time difference exceeds 5 seconds, since this is the max time.

ArcPy find points that are paired in time and space

3 Answers3