ArcPy Performance – Improving Performance with ArcGIS Cursors in Python

arcpycursorindexingperformance

I have a pretty big point feature class in a file geodatabase (~4 000 000 records). This is a regular grid of points with a 100m resolution.

I need to perform a kind of generalization on this layer. For this, I create a new grid where each point lies in the middle of 4 "old" points:

 *     *     *     *
    o     o     o
 *     *     *     *
    o     o     o
 *     *     *     *

[*] = point of the original grid – [o] = point of the new grid

The attribute value of each new point is calculated based on the weighted values of its 4 neighbors in the old grid. I thus loop on all the points of my new grid and, for each of them, I loop on all the points of my old grid, in order to find the neighbors (by comparing the values of X and Y in the attribute table). Once 4 neighbors have been found, we get out of the loop.

There is no methodological complexity here but my problem is that, based on my first tests, this script will last for weeks to complete…

Do you see any possibility to make it more efficient? A few ideas on the top of my head:

Index the fields X and Y => I did that but didn't notice any significant performance change
Do a spatial query to find the neighbors rather than an attribute-based one. Would that actually help? What spatial function in ArcGIS should do the job? I doubt that, e.g., buffering each new point will prove more efficient
Transform the feature class into a NumPy Array. Would that help? I haven't worked a lot with NumPy so far and I wouldn't like to dive into it unless someone tells me it might really help reducing the processing time
Anything else?

Best Answer

What if you fed the points into a numpy array and used a scipy cKDTree to look for neighbors. I process LiDAR point clouds with large numbers of points (> 20 million) in several MINUTES using this technique. There is documentation here for kdtree and here for numpy conversion. Basically, you read the x,y into an array, and iterate over each point in the array finding indices of points within a certain distance (neighborhood) of each point. You can use these indices to then calculate other attributes.

Related Solutions

ArcPy – Why Is Data Access Cursor Performance Enhanced in Recent Versions?

One of the developers of arcpy.da here. We got the performance where it is because performance was our primary concern: the main gripe with the old cursors were that they were slow, not that they lacked any particular functionality. The code uses the same underlying ArcObjects available in ArcGIS since 8.x (the CPython implementation of the search cursor, for example, looks a lot like code samples like this in its implementation except, you know, in C++ instead of C#).

The main two things we did to get the speedup are thus:

Eliminate layers of abstraction: the initial implementation of the Python cursor was based on the old Dispatch/COM based GPDispatch object, which enabled one to use the same API in any language that could consume COM Dispatch objects. This means that you had an API that was not particularly well-optimized for any single environment, but it also meant that there were a lot of layers of abstraction for the COM objects to advertise and resolve methods at runtime, for example. If you remember before ArcGIS 9.3, it was possible to write geoprocessing scripts using that same clunky interface many languages, even Perl and Ruby. The extra paperwork an object needs to do to handle the IDispatch stuff adds a lot of complexity and slowdown to function calls.
Make a tightly integrated, Python specific C++ library using Pythonic idioms and data structures: the idea of a Row object and the really strange while cursor.Next(): dance were just plain inefficient in Python. Fetching an item from a list is a very fast operation, and simplifies down to just a couple of CPython function calls (basically a __getitem__ call, heavily optimized on lists). Doing row.getValue("column") by comparison is more heavyweight: it does a __getattr__ to fetch the method (on which it needs to create a new bound method object), then call that method with the given arguments (__call__). Each part of the arcpy.da implementation is very closely integrated with the CPython API with a lot of hand-tuned C++ to make it fast, using native Python data structures (and numpy integration, too, for even more speed and memory efficiency).

You'll also notice that in nearly any benchmark (see these slides for example), arcobjects in .Net and C++ are still over twice as fast as arcpy.da in most tasks. Python code using arcpy.da is faster, but still not faster than a compiled, lower-level language.

TL;DR: da is faster because da is implemented in straight-up, unadulterated Arcobjects/C++/CPython which was specifically designed to result in fast Python code.

[GIS] How to improve performance of nested search cursors

I have to keep posting that I will never use embedded search cursors ever for relating tables to each other. Only dictionaries should be used to do what you are doing. Performance will increase 100 fold over a 3 level set of embedded search cursors, since dictionaries use random access from data stored in memory which is immediate, verses cursor's and sql that access the table linearly by repeated reading the data from disk which is incredibly sloooow. See my Blog on Turbo Charging Data Manipulation with Python Cursors and Dictionaries

Load each of the embedded cursors to their own dictionary using the related field as the key value using a style shown below that add all other fields into a list for the values. Then just run a cursor on the master data source and look up the key value in the first related table. Use a for loop on the first related table's value list to get the look up for the second level table key and then append the returned value of that table to the list you are building.

The conditional statements are virtually identical to what you have in your code and you should not think of dictionary logic any differently. The examples in my blog already show this kind of conditional logic. It occurs in the line that reads: if keyValue in valueDict:

Anyway, the logic of your code is less obvious than the logic of the dictionary code I use. I have tried to equate what you have done to how it should be done using dictionaries, but I may not have kept track of the list for loops that I needed to build. I made the Name Code1 value key of the two embedded cursors a tuple dictionary key for the most embedded level to make the matching easier. Not sure why you have origin and destination at that level, since they are not used. The second dictionary is basically using (Name, Code1) key as the equivalent of a multi-field composite key in a table relationship (something my blog discussed in general terms, but did not show in an example).

I have revised the code to use Paul's suggestion of using the defaultdict(list). The code has not been tested so some debugging is still most likely needed, but it illustrates the overall approach.

from time import strftime  
  
print "Start script: " + strftime("%Y-%m-%d %H:%M:%S")  
  
import arcpy  
  
sourceFC1 = parameters[2].valueAsText  
  
sourceFieldsList1 = ["Code","Name","Code1"] 
  
valueDict1 = defaultdict(list)  
with arcpy.da.SearchCursor(sourceFC1, sourceFieldsList1) as searchRows:  
    for searchRow in searchRows:  
        keyValue = searchRow[0]  
        valueDict1[keyValue].append((searchRow[1], searchRow[2]))  
                        
sourceFC2 = parameters[0].valueAsText  
  
sourceFieldsList2 = ["Name","Code1","Origin","Destination"] 
  
valueDict2 = defaultdict(list)  
with arcpy.da.SearchCursor(sourceFC2, sourceFieldsList2) as searchRows:  
    for searchRow in searchRows:  
        keyValue = (searchRow[0], searchRow[1])
        valueDict2[keyValue].append(searchRow[2:])
  
updateFC = parameters[1].valueAsText 
  
updateFieldsList = ["Station","Code"]  
  
with arcpy.da.UpdateCursor(updateFC, updateFieldsList, where_clause = sql_exp1) as updateRows:  
    for updateRow in updateRows:  
        if updateRow[0] != None:
        # store the Join value of the row being updated in a keyValue1 variable  
            keyValue1 = updateRow[1]  
            # verify that the keyValue1 is in the Dictionary  
            if keyValue1 in valueDict1:  
                # transfer the values stored under the keyValue1 from the dictionary to the updated fields.  
                for keyValue2 in valueDict1[keyValue1]: 
                    print(keyValue2)
                    if keyValue2 in valueDict2: 
                        listnumber.append(keyValue2[0])
                        listcodeFinal.append(keyValue2[1])
            codelinestr = ','.join(str(b) for b in linecodeFinal)
del valueDict1
del valueDict2  
  
print "Finished script: " + strftime("%Y-%m-%d %H:%M:%S")  
print "Finished script: " + strftime("%Y-%m-%d %H:%M:%S")

Best Answer

Related Solutions

ArcPy – Why Is Data Access Cursor Performance Enhanced in Recent Versions?

[GIS] How to improve performance of nested search cursors

Related Question