How to divide a feature class into ranges?

Question

I have a shapefile consisting of around 4000 polygons. The polygons are created from buffering around a similar number of islands somewhere in the Norwegian archipelago. Due to an error in the dissolve functionality using arcpy dissolve fails. I am assuming this is because arcpy.Dissolve_management can not handle too many overlapping features.

To try to work around this I plan on splitting my initial file into ten separate files consisting of around 400 features which are then dissolved separately. These ten files are then merged together and a final dissolve is performed.

How can I do such a selection? The arcpy.select_analysis allows for a where clause. How can I use this to select the n first features, then the n+1 to n*2 features, then the n*2+1 to n*3 features etc?

My final goal is to make a small arcpy function which could be a preliminary fix for us mere mortals relying on ESRIs bug fixing.

What exactly is the error you encounter? Have you considered simplifying the geometry? — blah238, Commented Nov 26, 2012 at 6:11
Simplified the geometry. Didnt help. It seems like a certain number of overlapping features leads to this error. Norwegian customer contact at geodata.no referred to it as: "NIM079373: Running a large number of features through the Dissolve or Buffer with dissolve option, hangs during process" Supposed to be fixed for 10.1 sp1, but they failed to do so. — ragnvald, Commented Nov 26, 2012 at 14:52
I have failed to find the referred to error code while googling it. Puzzled by this. It could be an ESRI internal reference. — ragnvald, Commented Nov 26, 2012 at 15:02

ccn · Accepted Answer · 2012-11-26 05:08:06Z

Use the looping and variable value as shown in Aragon's answer with the FID field (this is the zero-based object ID field for shapefiles) as the select field as shown in L_Holcombe's answer to generate the where clause, and all should be good. To program in the total number of features use the Get Count tool and divide by 10 assigned to the variable. Will also need Make Feature Layer. The lines with the '#' character need to be adjusted to suit your data and the '#' removed when ready to run.

#lyr = arcpy.MakeFeatureLayer_management(specify input file here)
#getcount = int(arcpy.GetCount_management(lyr).getOutput(0))
getcount = 4000     # just an example - remove line and use get count from above
numfiles = 10
numfeatures = getcount / numfiles

for i in range(0, getcount, numfeatures):
    select_features = specify output file here
    #dissolve_features = specify output file here
    j = i + numfeatures
    if i != 0:
        where_clause = '"FID" > ' + str(i) + ' AND "FID" <= ' + str(j)
    else:
        where_clause = '"FID" <= ' + str(j)
    print where_clause
    #arcpy.Select_analysis(lyr, select_features, where_clause)
    #arcpy.Dissolve_management(select_features, dissolve_features, {dissolve_field})

2 revs · Accepted Answer · 2012-11-26 17:36:22Z

I ended up making the function suggested below. In part inspired by inputs from other contributors on this question. The object handling is coarse and it creates a lot of temporary files.

I am sure this function could be made a lot better. If properly polished it could be a decent preliminary fix for the dissolve bug by ESRI. I am making it a community wiki.

# Import modules
import arcpy, os

def shape_dissolve(file_in, file_out,group_by):  
    count           =  0
    current_min     =  0
    current_max     =  group_by
    joinstring      =  []

    features_total    =  int(arcpy.GetCount_management(file_in).getOutput(0))

    while (current_max < (features_total+group_by)):

        resulting_file  = "%sbuffer_result_%s.shp" % (path_maps_basis,count)

        where_clause    = '"FID"> %s AND "FID" <= %s' %(current_min,current_max)

        arcpy.Select_analysis(file_in, resulting_file, where_clause)

        resulting_file_d = "%sbuffer_result_d_%s.shp" % (path_maps_basis,count)

        arcpy.Dissolve_management(resulting_file, resulting_file_d,"","","SINGLE_PART","")

        # delete temporary files
        arcpy.Delete_management(resulting_file)

        joinstring.append(resulting_file_d)

        current_min = group_by*count
        current_max = current_min + group_by

        count +=1

    resultbuffer = "%sresultbuffer.shp" % (path_maps_basis)

    arcpy.Merge_management(joinstring, resultbuffer,"")

    # delete temporary files
    for shapefile in joinstring:

        arcpy.Delete_management(shapefile)

    arcpy.Dissolve_management(resultbuffer, file_out,"","","SINGLE_PART","")

    # delete temporary files
    arcpy.Delete_management(resultbuffer)

#Calling the function with these sample values
group_by = 700
path_maps_basis     = "C:/"
multi_feature_file  = "%sbuffer_complex.shp"    % (path_maps_basis)
resulting_file      = "%sbuffer_dissolved.shp"  % (path_maps_basis)


shape_dissolve(multi_feature_file,resulting_file,group_by)

Just a note, this sequential grouping logic cannot be easily adapted to feature classes in a geodatabase because their OBJECTID values are not necessarily sequential (they can have gaps, and the highest OBJECTID can be greater than the feature count). The FID field in a shapefile is guaranteed to be sequential because it is generated on the fly and not stored in the shapefile: support.esri.com/es/knowledgebase/techarticles/detail/37480 Hence the slightly more complex partitioning logic used in my answer, but if you are only concerned with shapefiles then this should work fine. — blah238, Commented Nov 26, 2012 at 18:36
I am doing this one working with shapefiles, so it works for me. Not the ideal thing to do considering the vast options available in PostGIS. But quick and dirty does it this time :-) — ragnvald, Commented Nov 26, 2012 at 20:22

Community · Accepted Answer · 2017-05-23 12:39:42Z

4

You can use enumarete function in python.

myList = []

n = 99

for i, v in enumarete(myList):
    if i == n:
        print v 'first 100 item' #append your list
    if i > n and i < (n+2)*2:
        print v

Or you can use this script from SO:

def split_list(alist, wanted_parts=1):
    length = len(alist)
    return [ alist[i*length // wanted_parts: (i+1)*length // wanted_parts] 
             for i in range(wanted_parts) ]

A = [0,1,2,3,4,5,6,7,8,9]

print split_list(A, wanted_parts=1)
print split_list(A, wanted_parts=2)
print split_list(A, wanted_parts=8)

i hope it helps you...

edited May 23, 2017 at 12:39

CommunityBot

1

answered Nov 25, 2012 at 20:24

urcm

22.6k4 gold badges58 silver badges109 bronze badges

Thankyou - your code does the handling of the logic nicely. But the retrieval of the n features using arcpy still remains as a challenge.
– ragnvald
Commented Nov 25, 2012 at 20:35
rows = arcpy.UpdateCursor(full_buffer_temp) can be used to iterate through the features. So I need to push the features into a separate file. Close it and start anew.
– ragnvald
Commented Nov 25, 2012 at 20:40

Add a comment |

L_Holcombe · Accepted Answer · 2012-11-25 21:57:38Z

3

Perhaps something simple would work? Assuming "gid" is a column with consecutive integers:

arcpy.Select_analysis("in_features", "out_features", '"gid" < 400')
arcpy.Select_analysis("in_features", "out_features", '"gid" < 800 AND "gid" >= 400')

etc.

I haven't tried this, but it should work.

ArcGIS SQL expression reference: http://help.arcgis.com/en/arcgisdesktop/10.0/help/index.html#//00s500000033000000

answered Nov 25, 2012 at 21:57

L_Holcombe

1,3191 gold badge11 silver badges14 bronze badges

Add a comment |

Community · Accepted Answer · 2017-04-13 12:33:48Z

You can adapt the partitioning methodology used in this script for your purposes: Port “Producing Building Shadows” Avenue code to ArcGIS 10. Once you've got it working you could even experiment with multiprocessing by adapting other parts of the script (its purpose was to improve performance using multiple CPU cores).

The relevant function:

def getOidRanges(inputFC, oidfield, count):
    oidranges = []
    if procfeaturelimit > 0:
        message("Partitioning row ID ranges ...")
        rows = arcpy.SearchCursor(inputFC, "", "", oidfield, "%s A" % oidfield)
        minoid = -1
        maxoid = -1
        for r, row in enumerate(rows):
            interval = r % procfeaturelimit
            if minoid < 0 and (interval == 0 or r == count - 1):
                minoid = row.getValue(oidfield)
            if maxoid < 0 and (interval == procfeaturelimit - 1 or r == count - 1):
                maxoid = row.getValue(oidfield)
            if minoid >= 0 and maxoid >= 0:
                oidranges.append([minoid, maxoid])
                minoid = -1
                maxoid = -1
        del row, rows
    return oidranges

Note: procfeaturelimit is a global variable that specifies the maximum number of features each "partition" should be assigned.

Which is used like this:

count = int(arcpy.GetCount_management(inputFC).getOutput(0))
# Get the appropriately delmited field name for the OID field
oidfielddelimited = arcpy.AddFieldDelimiters(inputFC, oidfield)
oidranges = getOidRanges(inputFC, oidfield, count)
for o, oidrange in enumerate(oidranges):
    # Build a where clause for the given OID range
    whereclause = "%s >= %d AND %s <= %d" % (oidfielddelimited, oidrange[0], oidfielddelimited, oidrange[1])
    # Do your Select_analysis and whatever else here with the given whereclause

J Graham · Accepted Answer · 2015-01-20 04:52:38Z

Slight variation on @blah238's really great answer. Updated with arcpy.da goodness with appropriate tips'o'the cap documented:

The relevant function:

def getOIDChunks(inputFC, oidfield, groupsize):
    oids = []
    oidchunks = []

    #from https://arcpy.wordpress.com/2012/02/01/create-a-list-of-unique-field-values/
    with arcpy.da.SearchCursor(inputFC,['OID@']) as rows:
        oids = sorted({row[0] for row in rows})

    #originally from http://stackoverflow.com/a/312464/655787
    for i in xrange(0, len(oids), groupsize):
        oidchunks.append(oids[i:i+groupsize])

    return oidchunks

Use the function like this to get the where_clause. I did away with the count and just made the group size a user defined parameter like so:

# Get the appropriately delmited field name for the OID field
oidfielddelimited = arcpy.AddFieldDelimiters(inputFC, "OBJECTID")
oidranges = getOIDChunks(inputFC, "OBJECTID", 300)
for thisrange in oidranges:
    # Build a where clause for the given OID range
    whereclause = "%s >= %d AND %s <= %d" % (oidfielddelimited, min(thisrange), oidfielddelimited, max(thisrange))
    # Do your Select_analysis and whatever else here with the given whereclause
    print "currently working on chunk:  " + whereclause

Stack Exchange Network

How to divide a feature class into ranges?

6 Answers 6

Your Answer

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

Linked

Hot Network Questions

How to divide a feature class into ranges?

6 Answers 6

Your Answer

Sign up or log in

Post as a guest

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

Linked

Related

Hot Network Questions

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