1

I thought this would be a simple problem to solve, by my googling hasn't provided the answers so far.

I have created a near table (Using the 'Generate Near Table' in Arc) comparing two sets of points, and I have included the option to get the 10 closest points for each of the input points. So for each of the input points I have 10 near features(this can be see in the IN_FID column).

near table

I wish to calculate the average distance, plus the average angle (which I have converted to a bearing 0-360) of the 10 points for each unique value in IN_FID. I use the summary statistics tool with the case field set to IN_FID to get the average distance, however to get the average bearing I need to run a calculation that takes into account circular statistics.

My solution thus far within ModelBuilder, is to use the 'Row Selection' iterator with the IN_FID set as the 'Group by Field option', I then run the calculate field on the output for my average bearing. So for each unique value in the IN_FID column, I have one average bearing (from 10 points). However, this process is very slow, it takes about 3 seconds for each unique IN_FID, and with 10,000 points, is getting on for 8 hours of processing.

The summary statistics tool works in a matter of seconds, so I was wondering if there is a python script I can write and run in the calculate field tool code block, that means I don't have to use the 'Row Selection' iterator and speed up this process? (however, I'm not even sure this will speed it up!). My python isn't up to much at the moment, but my research has lead me to think I need to use search or update cursor, and some (nested?) for loops, but I don't have the knowledge to apply and write this yet, but I am learning using various examples I find online.

Does anyone have any python script examples that do a similar task that I could use as a basis for this task or, can point me in the direction of some useful tutorials?

Or do you have some alternative methods to overcoming this problem that I haven't thought of?

  • I think this would be done easiest and quickest using arcpy.da and a Python dictionary. ModelBuilder makes a poor choice for anything that is highly iterative - hard to configure and slow. The Python Parser of the Field Calculator is best suited to within row calculations. – PolyGeo Nov 18 '15 at 23:48
  • With average bearing it is tricky, because average of 1 and 359 = 180. Total cosine and sine use atan to define it – FelixIP Nov 19 '15 at 0:30
3

Never use nested cursor loops. Never ever! The concept of nested cursor loops is completely unnecessary.

Use a dictionary to collect information from a search cursor for each IN_FID as the key. The look up of the IN_FID from the dictionary is virtually instantaneous and should deal with 10,000 points in 1 minute or less.

Dictionaries can be populated from a search cursor to do statistical operations as each record is being read, or you can first collect all of the values associated with each FID key into a list within the dictionary and then use a for loop to do ordered statistical operations after sorting the list. These in memory processes are extremely fast. The summary value or sets of values can be stored back to the dictionary. Finally, you write all of your final statistical values using an update cursor that looks up each IN_FID from the dictionary to get the summary value you need to write, or you can create a new table and populate it using an insert cursor.

See my blog on Turbo Charging Data Manipulation Using Python Cursors and Dictionaries In particular look at the section on Using a Python Dictionary Built using a da SearchCursor to Replace a Summary Statistics Output Table.

I would need to know what circular statistical process you are doing with your distances and angles to come up with a modification to the script that would suit your exact needs. But I am sure that the basic outline of the script is what you need and can be modified to do what you want, and that this will meet your performance expectations.

  • Richard, many thanks for this, a great resource which helped me come up with a solution, which I have posted as an answer. However, the code still runs quite slow (0.5 seconds per unique IN_FID), which is still faster than my previous solution however!. Nevertheless, looking at the code, do you have any advice to improve the efficiency of the code? – meither Nov 20 '15 at 15:41
  • This code should be run as a stand alone script without a Field Calculation or a ModelBuilder model. Both are a waste of time. You can do everything faster in a script and eliminate the redundant use of the field calculator and the painfully slow model process. Rewrite the entire model as a script and use only cursors and dictionaries. A script can do everything the model does, probably in less than 1 minute. – Richard Fairhurst Nov 20 '15 at 19:06
  • The cursors are independent of the field calculator and are processing the entire point feature class everytime the calculation gets called. By running this code within the field calculator you are literally doing the average 10000 times for the entire 10,000 point feature class, which is 9,999 times too many. No wonder it is taking so long. To do all of your points takes only 5 seconds (the time it takes to do the calculation on just one selected point). Run the code just once in the model. – Richard Fairhurst Nov 20 '15 at 19:15
0

Using the post of Richard Fairhurt and this post for the bearing average I came up with the following solution:

import math
import arcpy
fc= "e:/NCCA/Temp/Temp.gdb/near_table"
bearings = ["IN_FID", "bearing", "bearing_avg"]
#creates a dictionary with the IN_FID as the key, and the bearing as the values
valueDict = {}
with arcpy.da.SearchCursor(fc, bearings) as searchRows:  
    for searchRow in searchRows:  
        keyValue = searchRow[0]  
        if not keyValue in valueDict:
            valueDict[keyValue] = [searchRow[1]] 
        else:  
            valueDict[keyValue].append(searchRow[1]) 

#calculates the bearing average of the values for each dictionary key and puts the average bearing into a new dictionary
averageBearingDict = {}
for key in valueDict:
    cosSum = 0.0
    sinSum = 0.0
    for bearingVal in valueDict[key]:
        theCos = math.cos(math.radians(float(bearingVal)))  
        theSin = math.sin(math.radians(float(bearingVal)))  
        cosSum += theCos  
        sinSum += theSin  
    N = len(valueDict[key]) 
    C = cosSum/N  
    S = sinSum/N  
    theMean = math.atan2(S,C)  
    if theMean < 0.0:  
        theMean += math.radians(360.0)  
    theMean_deg = math.degrees(theMean)

    averageBearingDict[key] = theMean_deg

#this puts the calculated bearing average back into the table - in row "bearing_avg"
with arcpy.da.UpdateCursor(fc, bearings) as rows:
    for row in rows:
        row[2] = averageBearingDict[row[0]]
        rows.updateRow(row)

This has been run as a stand alone script and completes in less than a minute.

  • I ran the code as a standalone script on a point feature class with 10000 points and it took only 10 seconds to run the entire script (7 second to read all of the 10000 point records into the dictionary, 1 second to do the averaging, and 2 seconds to write back to the 10000 points). I had to add an import arcpy line at the beginning and I removed the def calc() line and dedented the code below it. You should select all IN_FID records you want to process and only run the script once. Using this code with any other loop that selects IN_FID record groups first defeats the purpose of the code. – Richard Fairhurst Nov 20 '15 at 19:00
  • Richard - spot on again. Ran as a stand alone script, completes very quickly. Thanks very much for your help. My first proper python script works great! – meither Nov 22 '15 at 13:11

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.