Creating multipart polygons with ArcPy? Dissolve Failure?

6015
14
Jump to solution
08-24-2015 01:56 PM
ChadCordes
New Contributor II

I need the functionality of the editor > merge command for use in a python script tool. I realize this function is not directly accessible through ArcPy, but I am hoping someone can help with a work around?

The Issue: I have a series of single part polygon features which do not touch or overlap, but have common field values. I need all these features to exist as multipart features for the random point generation to operate properly and consider all these regions equally. However, the dissolve function inconsistently fails leaving single part polygons.

Essentially, I am trying to go from single part geometries to multipart by creating feature layer selections and dissolving on that feature layer selection.

I have included a code snippet below. The code works and does exactly as required, but it is inconsistent. The dissolve fails when preceded by other scripts or when significant edits with debugging have been made to other parts of this script. However, when called from a new .mxd document it runs consistently well.

Any advice and support would be gratefully accepted!
Thanks!

#Custom exceptions.
class LicenseException(Exception):
    pass


#Import system modules.
import arcpy
import sys
import os
import traceback
from arcpy import env
from arcpy.sa import *


#Enable geoprocessor logging.
arcpy.SetLogHistory(True)


try:
    #---Check the installed ArcGIS version.---
    iv = arcpy.gp.GetInstallInfo()


    version = " ArcGIS %s : %s" % ('Version', iv.get('Version','Version Unknown'))


    if version.find("10.") > 0:
        ArcGIS10 = True
    else:
        ArcGIS10 = False
        arcpy.AddError("ERROR: This tool requires ArcGIS version 10.1 or greater...exiting script.")
        sys.exit("")


    del iv, version, ArcGIS10


    #---Check the ArcGIS administrative license.---
    if arcpy.CheckProduct("ArcInfo") == "AlreadyInitialized" or arcpy.CheckProduct("ArcInfo") == "Available":  
        pass
    else:
        raise LicenseException()


    #---Check out extension licenses.---
    if arcpy.gp.CheckExtension("spatial") == "Available":
        arcpy.CheckOutExtension("spatial")
    else:
        arcpy.AddError("ERROR: Spatial Analyst extension is unavailable. Please enable Spatial analyst and try again...exiting script.")
        sys.exit("")


    #---Get tool share folder file path.---
    mxd = arcpy.mapping.MapDocument("Current")
    mxdpath = mxd.filePath
    (filepath, filename) = os.path.split(mxdpath)
    MCP_Path = str(filepath)


    #---Set the geoprocessing environment.---
    workspace = MCP_Path + "\OutputData\ToolOut.gdb"
    scratchspace = MCP_Path + "\OutputData\ScratchFolder\Scratch.gdb"


    env.scratchWorkspace = scratchspace
    env.workspace = workspace
    env.overwriteOutput = True


    #---Read user input.---
    BI = arcpy.GetParameterAsText(0)


    BIdesc = arcpy.Describe(BI) 


    #---Execute Reclass.---
    DS_RC_filename = BIdesc.name + "_DS_RC"
    DS_RC_out = os.path.join(workspace, DS_RC_filename)


    DS_Reclass = Reclassify(BI, "Value", 
    RemapRange([[-1000,-75,8],[-75,-50,7],[-50,-35,6],[-35,-20,5],[-20,-12,4],[-12,-6,3],[-6,-3,2],[-3,-1,1],[-1,0,"NODATA"]]), "NODATA")
    DS_Reclass.save(DS_RC_out)


    #---Execute raster to polygon conversion.---
    DS_Polygon_filename = BIdesc.name + "_DS_Polygon"
    DS_Polygon = os.path.join(workspace, DS_Polygon_filename)


    arcpy.RasterToPolygon_conversion(DS_RC_out, DS_Polygon, "SIMPLIFY", "VALUE")


    if arcpy.Exists(DS_RC_out):
        arcpy.Delete_management(DS_RC_out)
    del DS_RC_filename, DS_Polygon_filename 


    #---Select polygons of areas >= 1000m2.---
    DS_PS_out_filename = BIdesc.name + "_DS_PSelect"
    DS_PS_out = os.path.join(workspace, DS_PS_out_filename)


    sql = '{0} >= {1}'.format("Shape_Area",'1000')
    if arcpy.Exists("DS_PS"):
        arcpy.Delete_management("DS_PS")
    DS_PS = arcpy.management.MakeFeatureLayer(DS_Polygon,"DS_PS",sql)
    arcpy.CopyFeatures_management("DS_PS", DS_PS_out) #Done to hopefully stabilize the following selection


    arcpy.Delete_management(DS_PS)
    del DS_PS_out_filename, sql


    #---------------------Perform Dissolve----------------------------------    
    DS_F_out_filename = BIdesc.name + "_DS"
    DS_F_out = os.path.join(workspace, DS_F_out_filename)


    arcpy.Dissolve_management(DS_PS_out,DS_F_out,"grid_code","#","MULTI_PART","DISSOLVE_LINES") 


    del BsMmaxS, sql


    #---Check to ensure dissolve function worked properly.---
    sql = '{0} = {1}'.format("grid_code",1)
    arcpy.MakeTableView_management(DS_F_out, "DS_F_out_TV",sql)
    count = int(arcpy.GetCount_management("DS_F_out_TV").getOutput(0))
    if count > 1:
        del sql
        if arcpy.Exists(DS_F_out):
            arcpy.Delete_management(DS_F_out)
        arcpy.Delete_management("DS_F_out_TV")
        arcpy.AddError("Geoprocessing Error due to memory limitations with 32bit processing. Assuming 64bit processing is unavaiable, please close your ArcMap instance, open it again and reprocess...exiting script.")
        #Above error prompt is my current best guess since the polygons can be quite complex and can contain many verticies. 
        sys.exit("")
    else:
        arcpy.Delete_management("DS_F_out_TV")
        del sql


    #Continue code here..
Tags (3)
0 Kudos
1 Solution

Accepted Solutions
curtvprice
MVP Esteemed Contributor

I am seeing a similar deal with Dissolve run at the end of a large script. The dissolve works fine from ArcMap and run from a clean Python session with similar environment. Here is my error message:

21:18:53 43.63 Dissolve
ERROR: File "Q:\tools\nhr_raster\HRStepL.py", line 1160, in RasterProc
    "OBJECTID COUNT", "MULTI_PART")
Executing: Dissolve lyrCat C:\Working\StepL_0604_20151123\work.gdb\vpucat_poly1 GRIDCODE "OBJECTID COUNT" MULTI_PART DISSOLVE_LINES
Start Time: Mon Nov 23 21:18:56 2015
Sorting Attributes...
Dissolving...
ERROR 999999: Error executing function.
ERROR: Invalid Topology [Cannot open cache file.]
ERROR: Failed to execute (Dissolve).
Failed at Tue Nov 24 01:23:57 2015 (Elapsed Time: 4 hours 5 minutes 1 seconds)

Here's what happens when run from a straight Python command line:

Executing: Dissolve C:\Working\StepL_0604_20151123\work.gdb\vpucat_poly0 C:\Working\StepL_0604_20151123\work.gdb\vpucat_poly0a Grid
ode "objectid count" MULTI_PART DISSOLVE_LINES
Start Time: Tue Nov 24 10:27:58 2015
Sorting Attributes...
Dissolving...
Succeeded at Tue Nov 24 10:30:09 2015 (Elapsed Time: 2 minutes 11 seconds)

I have set ARCTMPDIR to a writeable folder (Esri KB 25995) and have even tried running half the polygons and get the same behavior. Given how things are hanging before they crash i suspect there is a resource issue here (ie lack of RAM) and that things eventually time out.

I don't get the problem (I don't think) when running 64 bit, which is even more evidence that this is a resources issue. (Along with your experience that restarting the mxd lets it work - restarting the mxd would clear layers and other cruft from memory.)

I am going to go back and try deleting some layers from memory in the hope that will solve my problem. My scripts clean up memory in a finally block but I think there is enough stuff in this script that I need to do some pre-cleanup!

I'm also going to try to add some garbage collection to hopefully release some memory that Python may be holding onto. Garbage collection happens automatically but a complicated script may benefit from a gc.collect() before the max possible memory may be needed.

import gc
gc.collect()

UPDATE: I've started an incident with Esri support and they pretty much told me that Dissolve is not stable with really large inputs (my input that was failing is a Raster To Polygon output of Watershed() output with about 70,000 polygons). The workaround is to tile the data and dissolve smaller datasets, or run in x64 python (for example, x64 background processing). (Not an option for me, unfortunately.)

View solution in original post

14 Replies
curtvprice
MVP Esteemed Contributor

Are the offending polygons that won't dissolve godzillas?

Please format your code.

0 Kudos
ChadCordes
New Contributor II

Hi Curtis;

Great suggestions. However, I just checked the polygons and they're not godzillas. I just got a failure from a dissolve with vertex counts ranging from 30 to 120.

0 Kudos
JoshuaBixby
MVP Esteemed Contributor

What is the specific error message, can you post it please?

0 Kudos
ChadCordes
New Contributor II

There is no error message. The dissolve completes, but the multiple single part features are not dissolved into a mulitpart feature class.

0 Kudos
DavidWasserman
Occasional Contributor III

Hi Chad,

I actually had to deal with something similar. I guess my main question is what progress did you make either the  eliminate tool? It sounds like that is what you want to do in the sense it might allow you to query/exclude specific fields. Is the area issue the concern?

David

David Wasserman, AICP
0 Kudos
ChadCordes
New Contributor II

No the area issue is not the concern. The issue is I am trying to create multipart polygons from a series of single part features using ArcPy. The operation is similar to that performed with the editor > merge functionality. This however is not supported with ArcPy. I have found a work around by using feature layer selections and then dissolving on a given field, but the dissolve function is inconsistently failing for some unknown reason. It primarily fails when it has been preceded by other scripts, or when significant modifications have been made to other lines of code within the script which calls the dissolve function.

I appreciate the suggestion with the eliminate tool, but I don't think that will solve my problem. Thanks though.

0 Kudos
DavidWasserman
Occasional Contributor III

I misunderstood then. Thanks for the clarification.

David Wasserman, AICP
0 Kudos
JoshuaBixby
MVP Esteemed Contributor

From a purely Python programming perspective, having 100+ lines of code within a single try block is unusual.  What does your except clause look like?

0 Kudos
ChadCordes
New Contributor II

Granted it is a large script.

#---Exception handling.---
except SystemExit:
  pass

except KeyboardInterrupt:
  arcpy.AddMessage("Interruption requested...exiting script")

except arcpy.ExecuteError:
  msgs = arcpy.GetMessages(2)
  arcpy.AddError(msgs)

except LicenseException:
  arcpy.AddError("ERROR: The BsM 2015 Map Creation Package requires the ArcGIS Advanced(ArcInfo) license. Please see your ArcGIS License Administrator...exiting script.")

except:
  tb = sys.exc_info()[2]
  tbinfo = traceback.format_tb(tb)[0]
  pymsg = "PYTHON ERRORS:\nTraceback info:\n" + tbinfo + "\nError Info:\n" + str(sys.exc_info()[1])
  msgs = "ArcPy ERRORS:\n" + arcpy.GetMessages(2) + "\n"
  arcpy.AddError(msgs)
  arcpy.AddError(pymsg)

  if arcpy.Exists(DS_F_out):
       arcpy.Delete_management(DS_F_out)

finally:
  #Clean.
  if arcpy.Exists(DS_Polygon):
       arcpy.Delete_management(DS_Polygon)
  if arcpy.Exists(DS_PS_out):
       arcpy.Delete_management(DS_PS_out)

  #Check In Extensions.
  arcpy.CheckInExtension("spatial")

  #Refresh the Table of Contents and Active View windows.
  arcpy.RefreshTOC()
  arcpy.RefreshActiveView()

  #---Include text to improve batch process result readability.---
  arcpy.AddMessage("\n")
0 Kudos