Tuesday, February 23, 2016

ArcGIS Erase without an ArcGIS advanced license

After creating code for converting Polygons to Lines and Vertices to Points, I decided to try my hand at creating an Erase function that didn't require ArcGIS's advance license. The code I came up with can erase polygons from lines, lines from lines, and polygons from polygons. The inputs are:


  • inFc: The input feature class
  • eraseFc: The feature class whose features will be erased from the input feature class
  • outLoc: The output file geodatabase
  • outName: The final name of the output feature class
  • workspace: A file geodatabase workspace
  • notices: Set to True for printed progress information


The code:

def Erase (inFc, eraseFc, outLoc = "", outName = "", workspace = "in_memory", notices = False):
    import os
    addToToc = arcpy.env.addOutputsToMap
    arcpy.env.addOutputsToMap = False
    if not outName:
        outName = os.path.basename (inFc) + "_Erase"
    if not outLoc:
        catPath = arcpy.Describe (inFc).catalogPath
        outLoc = os.path.dirname (catPath)
    outFc = UniqueFileName (outLoc, outName)
    outLoc, outName = os.path.split (outFc)

    garbage = []

    inType = arcpy.Describe (inFc).shapeType
    eraseType = arcpy.Describe (eraseFc).shapeType

    #lines and lines
    if inType == "Polyline" and eraseType == "Polyline":
        #buffer both lines: flat, left
        inBuffFc = UniqueFileName (workspace)
        arcpy.Buffer_analysis (inFc, inBuffFc, ".01 FEET")
        garbage += [inBuffFc]
        eraseBuffFc = UniqueFileName (workspace)
        arcpy.Buffer_analysis (eraseFc, eraseBuffFc, ".01 FEET")
        garbage += [eraseBuffFc]
        
        #union input buffer with erase buffer
        unionFc = UniqueFileName (workspace)
        arcpy.Union_analysis ([inBuffFc, eraseBuffFc], unionFc)
        garbage += [unionFc]
        #select union output by erase buffer
        unionLyr = UniqueLayerName (unionFc)
        garbage += [unionLyr]
        arcpy.SelectLayerByLocation_management (unionLyr, "WITHIN", eraseBuffFc)
        #switch union selection
        arcpy.SelectLayerByAttribute_management (unionLyr, "SWITCH_SELECTION")
        #check for selection
        if not arcpy.Describe (unionLyr).FIDSet:
            arcpy.CreateFeatureclass_management (outLoc, outName, "POLYLINE", inFc, "", "", inFc)
            for trash in garbage:
                arcpy.Delete_management (trash)
            return outFc
        #clip line
        arcpy.env.addOutputsToMap = addToToc
        arcpy.Clip_analysis (inFc, unionLyr, outFc)
        

    #lines and polygons
    elif inType == "Polyline" and eraseType == "Polygon":
        #buffer both lines: flat, left
        inBuffFc = UniqueFileName (workspace)
        arcpy.Buffer_analysis (inFc, inBuffFc, ".01 FEET")
        garbage += [inBuffFc]
        
        #union input buffer with erase buffer
        unionFc = UniqueFileName (workspace)
        arcpy.Union_analysis ([inBuffFc, eraseFc], unionFc)
        garbage += [unionFc]
        #select union output by erase buffer
        unionLyr = UniqueLayerName (unionFc)
        garbage += [unionLyr]
        arcpy.SelectLayerByLocation_management (unionLyr, "WITHIN", eraseFc)
        #switch union selection
        arcpy.SelectLayerByAttribute_management (unionLyr, "SWITCH_SELECTION")
        #check for selection
        if not arcpy.Describe (unionLyr).FIDSet:
            arcpy.CreateFeatureclass_management (outLoc, outName, "POLYLINE", inFc, "", "", inFc)
            for trash in garbage:
                arcpy.Delete_management (trash)
            return outFc
        #clip line
        arcpy.env.addOutputsToMap = addToToc
        arcpy.Clip_analysis (inFc, unionLyr, outFc)

    #polygons and polygons
    elif inType == "Polygon" and eraseType == "Polygon":
        unionFc = UniqueFileName (outLoc, outName)
        arcpy.env.addOutputsToMap = addToToc
        arcpy.Union_analysis ([inFc, eraseFc], unionFc, "NO_FID")
        arcpy.env.addOutputsToMap = False
        unionLyr = UniqueLayerName (unionFc)
        garbage += [unionLyr]
        arcpy.SelectLayerByLocation_management (unionLyr, "WITHIN", eraseFc)
        arcpy.DeleteFeatures_management (unionLyr)

        inFlds = [f.name for f in arcpy.ListFields (inFc)]
        delFlds = [f.name for f in arcpy.ListFields (unionLyr) if not f.name in inFlds and not
                   f.type in ("Geometry", "OID") and f.required == False]
        arcpy.DeleteField_management (unionLyr, delFlds)
        arcpy.env.addOutputsToMap = addToToc

    for trash in garbage:
        arcpy.Delete_management (trash)    
    
    return outFc

def UniqueFileName(Location = "in_memory", Name = "file", Extension = ""):
    """
    returns an unused file name
    'Location' will be the file folder
    'Name' will be the file name plus numeric extension if needed
    'Extension' is the file extension

    similar to arcpy.CreateUniqueName
    """
    import arcpy
    import os

    if Extension:

        outName = os.path.join (Location, Name + "." + Extension)
    else:
        outName = os.path.join (Location, Name)
    i = 0
    while arcpy.Exists (outName):
        i += 1
        if Extension:
            outName = os.path.join (Location, Name + "_" + str(i) + "." + Extension)
        else:
            outName = os.path.join (Location, Name + "_" + str(i))
        
    return outName

def UniqueLayerName(inFC = None, sql = None):
    """
    returns an unused feature layer name
    """
    import arcpy
    import os
    lyrName = os.path.join ("lyr0")
    i = 0
    while arcpy.Exists (lyrName):
        i += 1
        lyrName = "lyr" + str(i)
    if inFC:
        arcpy.MakeFeatureLayer_management (inFC, lyrName, sql) 
    return lyrName
profile for Emil Brundage at Geographic Information Systems Stack Exchange, Q&A for cartographers, geographers and GIS professionals

Wednesday, February 3, 2016

ArcGIS polygon to lines and feature vertices to points without an ArcGIS advanced license

I've started taking a little time to create python functions, using standard license code, that replicate some of ArcGIS's tools that are only available with an advanced license. So far, I've created functions for the Polygon To Line and the Feature Vertices To Points tools. Each make use of cursors and geometries to get the job done. The two functions take an input feature class's full catalog path, a workspace (file geodatabase or in_memory), and the file name for the new feature class, and return the full path to the function's output, which has been generated in the indicated workspace. Both also make use of my UniqueFileName function mentioned in a previous blog.

Which tools should I code next?

The code:

#polygons to lines
def PolygonToLines (inFc, workspace, fcName):

    outFc = UniqueFileName (workspace, fcName)
    try:
        #let esri do it (advanced license)
        arcpy.FeatureToLine_management (inFc, outFc)
        return outFc
    except:
        pass
    #no advanced license
    array = arcpy.Array ()
    sr = arcpy.Describe (inFc).spatialReference
    outPath, outName = os.path.split (outFc)
    arcpy.CreateFeatureclass_management (outPath,
                                         outName,
                                         "POLYLINE",
                                         spatial_reference = sr)
    with arcpy.da.InsertCursor (outFc, "SHAPE@") as iCurs:
        with arcpy.da.SearchCursor (inFc, "SHAPE@") as sCurs:
            for geom, in sCurs:
                for i in range (geom.partCount):
                    lastPart = None
                    part = geom.getPart (i)
                    for pnt in part:
                        if not pnt:
                            lastPart = None
                            continue
                        x = pnt.X
                        y = pnt.Y
                        if not lastPart:
                            lastPart = (x, y)
                            continue
                        thisPart = (x, y)
                        lastPoint = arcpy.Point ()
                        lastPoint.X = lastPart [0]
                        lastPoint.Y = lastPart [1]
                        thisPoint = arcpy.Point ()
                        thisPoint.X = thisPart [0]
                        thisPoint.Y = thisPart [1]
                        array.add (lastPoint)
                        array.add (thisPoint)
                        polyline = arcpy.Polyline (array)
                        array.removeAll ()
                        lastPart = thisPart
                        row = (polyline,)
                        iCurs.insertRow (row)
    del iCurs
    del sCurs
    return outFc


#feature vertices to points
def VerticesToPoints (inFc, workspace, fcName):

    outFc = UniqueFileName(workspace, fcName)
    try:
        #let esri do it (advanced license)
        arcpy.FeatureVerticesToPoints_management (inFc, outFc)
        return outFc
    except:
        pass
    #no advanced license
    sr = arcpy.Describe (inFc).spatialReference
    outPath, outName = os.path.split (outFc)
    arcpy.CreateFeatureclass_management (outPath,
                                         outName,
                                         "POINT",
                                         spatial_reference = sr)
    with arcpy.da.InsertCursor (outFc, "SHAPE@") as iCurs:
        with arcpy.da.SearchCursor (inFc, "SHAPE@") as sCurs:
            for geom, in sCurs:
                for i in range (geom.partCount):
                    part = geom.getPart (i)
                    for pnt in part:
                        if not pnt:
                            continue
                        row = (pnt,)
                        iCurs.insertRow (row)
    del iCurs
    del sCurs
    return outFc


#generate a unique ArcGIS file name
def UniqueFileName(location = "in_memory",
                   name = "file",
                   extension = ""):
    
    if extension:
        outName = os.path.join (location, name + "." + extension)
    else:
        outName = os.path.join (location, name)
    i = 0
    while arcpy.Exists (outName):
        i += 1
        if extension:
            outName = os.path.join (location,
                                    "{0}_{1}.{2}".format (name, i,
                                                          extension))
        else:
            outName = os.path.join (location,
                                    "{0}_{1}".format (name, i))
    return outName
profile for Emil Brundage at Geographic Information Systems Stack Exchange, Q&A for cartographers, geographers and GIS professionals