1

I have a numpy data processing loop that repeats ~50,000x and outputs ~23 data values. I want to write this data to the attribute table of my ~50,000-point multipoint feature class. With each loop, I generate 8 individual values, plus 5 sets of three values. (I'm still deciding whether to group them as 3 sets of 5 or 5 sets of 3, TBD by the answer to below.)

With this number of iterations, I want to be mindful of processing time. How should I approach organizing & writing the data? Should I:

  1. ...write 23 fields of data to the output file row-by-row at the end of each loop? (Or is the processing time of writing to file too long to do ~50Kx?)
  2. ...or collect them into 23 separate ~50K-length numpy arrays and write them column-by-column after all looping is done? (Or would I likely run into memory issues?)

Here's my approach for the second method:

#23 data columns to be added
fieldlist = ["azimuth", "orient", "crit#f", "freq1", "amp1", "pwr%1", "freq2", "amp2", "pwr%2", "freq3", "amp3", "pwr%3", "freq4", "amp4", "pwr%4", "freq5", "amp5", "pwr%5", "Nspurs", "len_mean", "len_sd", "ht_mean", "ht_sd"]

for fieldname in fieldlist:
    arcpy.AddField_management(gridcopy, fieldname, "FLOAT") #Add new fields
    with arcpy.da.UpdateCursor(gridcopy, [fieldname]) as cursor: #Add data
        i = 0
        for row in cursor:
            row = X[i]
            cursor.updateRow(row)
            i += 1
  • The slow part is your updateRow(), add all the fields and then calc each row. – Michael Stimson Mar 24 '17 at 0:02
  • @MichaelMiles-Stimson: So that means you're recommending option 1? If I follow, you're saying with option 1 I have to do the slow part ~50K times; but with option 2 I'll do it 23*50K times. Correct? – joechoj Mar 24 '17 at 3:08
  • 1
    Yes, update the data row by row not field in fields. Add all the fields to your database first (if they don't already exist) then calculate all the fields for each row and update. If you're really worried about how slow it is then use ArcObjects in C# (or even better C++) instead of python. – Michael Stimson Mar 24 '17 at 3:51
4

If you don't have any data already preserved in the output feature class, it might be faster to insert new rows instead of updating existing ones. For this, you would use da.InsertCursor.

If you still would like to update the rows, then do the following:

#23 data columns to be added
fieldlist = ["azimuth", "orient", "crit#f", "freq1", "amp1", "pwr%1", "freq2", "amp2", "pwr%2", "freq3", "amp3", "pwr%3", "freq4", "amp4", "pwr%4", "freq5", "amp5", "pwr%5", "Nspurs", "len_mean", "len_sd", "ht_mean", "ht_sd"]

for fieldname in fieldlist:
    arcpy.AddField_management(gridcopy, fieldname, "FLOAT") #Add new fields

with arcpy.da.UpdateCursor(gridcopy, fieldlist) as cursor: #Add data
    for row in cursor:
        row = X 
        cursor.updateRow(row)

~50K rows is nothing, will be updated in no time at all. Since you are working with numpy arrays, it might also be worth to check arcpy.da methods for converting numpy arrays into feature classes such as NumPyArrayToFeatureClass which will convert a NumPy structured array to a point feature class. They are very fast, I've been using arcpy.da module functions for large numpy arrays (millions of records) and they are very fast.

Performance information update:

It will be faster to use the UpdateCursor than to create geodatabase tables first from numpy arrays as creating new objects in a geodatabase does have an overhead and using ArcGIS-based joins would be very slow. Look at the execution time results - updating actual rows takes just 2 seconds; the slowest part is importing arcpy and getting count of features:

from functools import wraps
import numpy as np
import time

fields_list = ['field{}'.format(i) for i in xrange(1,24)]
fc = r'C:\GIS\Temp\ArcGISHomeFolder\Default.gdb\points'

#----------------------------------------------------------------------
def report_time(func):
    '''Decorator reporting the execution time'''
    @wraps(func)
    def wrapper(*args, **kwargs):
        start = time.time()
        result = func(*args, **kwargs)
        end = time.time()
        print(func.__name__, round(end-start, 3))
        return result
    return wrapper

@report_time
#----------------------------------------------------------------------
def setup():
    """set up sample data and add fields to feature class"""
    import arcpy
    total_features = int(arcpy.GetCount_management(fc).getOutput(0))
    print total_features
    #for f in fields_list:
        #arcpy.AddField_management(in_table=fc, field_name=f, field_type='FLOAT')
    data = np.random.uniform(low=0.5, high=13.3, size=(total_features, 23))
    return data

@report_time
#----------------------------------------------------------------------
def update_cursor(data):
    """update feature class fields using da.UpdateCursor"""
    with arcpy.da.UpdateCursor(fc, fields_list + ['OID']) as cursor: #Add data
        for row in cursor:
            row[0:23] = data[row[-1] - 1]
            cursor.updateRow(row)

#----------------------------------------------------------------------
def main():
    """"""
    data = setup()
    update_cursor(data)

if __name__ == '__main__':
    main()

49696 ('setup', 20.678) ('update_cursor', 1.837)

  • I do have data in the output feature class, but it's just the OBJECTID and Shape data for the 50K points in the multipoint file. So if I use AddFields to this file, I'll end up with a large table of empty cells that I'll have to backfill with data. I wonder if I can construct a separate table, and then perform a Join on the two? Or should I just go with the UpdateCursor? – joechoj Mar 24 '17 at 23:37
  • 1
    Please see the update in the answer body – Alex Tereshenkov Mar 26 '17 at 10:46
  • Thanks very much for this demo. I'm not sure what's in fc: does this demo go through ~50K rows? I really like your suggestion of NumpyArrayToFeatureClass: that actually seems perfect. Thanks for drawing my attention to that. Also, I can't figure out what's going on in lines 2&4 of update_cursor. In line 4, you're assigning values to the row's fields. row[-1] is the last field in the row; why subtract 1 from that? Line 2: Why add ['OID'] to fields_list? – joechoj Mar 30 '17 at 3:40
  • 1
    No problem at all! Yes, fc is actually a fc with 50K rows. The fc contains those fields in the fields_list but also OID field. The numpy array contains the values for 23, but not OID value. However, the array has an order - so the first item in the array (Python index of 0) should be populating fields for feature with OID 1. So, data[row[-1] - 1] would give me the first item in the numpy array, because row[-1] would give me OID value of the feature. OID starts with 1 in ArcGIS, however Python index is 0-based. Thus I have to subtract 1. – Alex Tereshenkov Mar 30 '17 at 6:34
  • I'm still having trouble getting your suggested solution to work. I'm curious why you eliminated the 'i' variable I had included in my example. In your 2nd-last line you had row = X, but don't you need a means of controlling which row in 'X' gets written to 'row'? – joechoj May 4 '17 at 22:02

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.