How to improve performance when using ArcGIS cursors in Python with big tables?

Question

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?

Answer 1

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.

Answer 2

I am with Barbarossa... arcpy cursors are insanely lame, so I only use them to traverse a table or feature class exactly one time. If I can't get the job done in one cycle, I use the cursor to fill up some other kind of data structure and work with that.

If you do not want to hassle with numpy just make a simple python dictionary where you use your coordinates as a simple text key, and fill in the attributes you need for calculation into a list as the value of the dictionary item.

In a second step you can easily get the values you need to calculate a point by simply getting them from your dictionary (which is incredibly fast, because of the dictionaries hashindex of items).

Answer 3

For a regular grid, it should be by far more efficient to work in a raster format. Convert your first grid into a raster, the you can resample at the same resolution using a bilinear interpolator but shifting your output image by 1/2 pixel in X and Y, and back again to points if you still need to have points.

EDIT : for complex decisions rules, you can convert each of the fields that you need as a new raster band, then you make four copies of those bands and you shift you raster in the 4 directions by 1/2 pixel (+50, -50), (+50,+50), (-50,-50) and (-50,+50). Then you can use regular map algebra

Answer 4

Thanks everybody for your help!

I finally found a very non-pythonic way to solve this issue... What was actually taking the most computing time was to find the 4 neighbors of each point. Rather than using the X and Y attributes (either with an arcpy cursor or within another data structure, such as a python ditionary), I ended up using the ArcGIS tool Generate near table. I assume this takes advantage of the spatial indexes and the performances are obviously much much higher, without me having to implement the index myself.

Answer 5

Problem with cursors is that you can cycle through them in one way only and you cant go back. Although not recommended, you can populate the feautres into a structure if you are planning to revisiting them.

If you were able to process your features in a single loop, I suggest enabling recycling. It is a parameter on your search featureclass function that lets python to reuse the memory allocated by old features and make traversing the features in a cursor much faster. You can process your grid 80% faster.

Problem is that you can't enable recycling if you are planning to store retrieved features from a cursor.