Python Code issue

Question

Original User: luke.kaim

Hi Everyone,

I am trying to write python code to detect if there is overlapping polygon from two feature classes. The goal of this code is to get overlapping polygons and there area. I start by trying to get the area of the polygons themselves. This I am trying to do using arcpy.geometry. Next I want to use the overlap geometry operation in 10.1. This is a boolean operation that will speed up the intersection command later. Then I do the intersection and the symmetric difference. The end goal is to get polygon overlap/difference.

I think I have a few bugs in my code, but I am struggling to get the logic and the syntax correct. I tried to model my code off of this thread: http://forums.arcgis.com/threads/60968-area-of-overlap-between-polygons-in-two-different-feature-classes?highlight=arcpy.Geometry%28%29 .

Any help would be greatly appreciated.
Code:def main():
    try:
        import arcpy, sys, traceback, os, string, locale
        arcpy.env.overwriteOutput = True

from arcpy import env
        env.workspace =r"C:\Users\Luke Kaim\Documents\University of Maine\Fall_2012\Volunteer Geographic Information\simalar"
        
        VGIFile = r"C:\Users\Luke Kaim\Documents\University of Maine\Fall_2012\Volunteer Geographic Information\simalar\polygon.shp"
        VGIFileShapeName = arcpy.Describe(VGIFile).ShapeFieldName
        VGIRows = arcpy.SearchCursor(VGIFile)

GNIS_FC = r"C:\Users\Luke Kaim\Documents\University of Maine\Fall_2012\Volunteer Geographic Information\simalar\GNIS.shp"
        GNIS_FCShapeName = arcpy.Describe(GNIS_FC).ShapeFieldName
        GNIS_FCRows = arcpy.SearchCursor(GNIS_FC)

#Create an empty Geometry object
        g = arcpy.Geometry()
        area=0
        # get Geometry of every VGI polygon. Then get the area of each polygon
            for  VGIFileShapeName in VGIRows:
                feat =  VGIFileShapeName.getValue(shapeName)
                area += feat.area
print area

#Create an empty Geometry object
        g = arcpy.Geometry()
        area=0
        # get Geometry of every GNIS polygon. Then get the area of each polygon
             for  GNIS_FCShapeName in GNIS_FCRows:
                feat =  GNIS_FCShapeName.getValue(shapeName)
                area += feat.area
print area

for VGIRow in VGIRows:
                    VGIGEOM = VGIRow.getValue(VGIFileShapeName)

GNIS_Rows = arcpy.SearchCursor(GNIS_FC)
                    for GNIS_Row in GNIS_Rows
                        GNIS_Geom = GNIS_Row.getValue(GNIS_FCShapeName)
                        if GNIS_Geom.overlaps(VGIGEOM):
                            arcpy.Intersect_analysis([GNIS_Geom,VGIGEOM],"in_memory/Intersection")
                            arcpy.SymDiff_analysis([GNIS_Geom,VGIGEOM],"in_memory/SymDiff")
                            
                            intersectionRows = arcpy.SearchCursor("in_memory/Intersection")
                    for intersectionRow in intersectionRows:
                        intersectionArea = intersectionRow.Shape.area
                        print "Overlap Area: " + str(intersectionArea)

intersectionRows2 = arcpy.SearchCursor("in_memory/"in_memory/SymDiff"")
                    for intersectionRow1 in intersectionRows2:
                            intersectionArea2 = intersectionRow1.Shape.area
                        print "Overlap Area :in_memory/SymDiff " + str(intersectionArea2)

del   VGIRows, GNIS_FCRows,GNIS_Rows

Thank you,
Luke Kaim

Thank you
Luke Kaim (SIE)
Lucas.Kaim@maine.edu
(914)263-7866
"Don�??t complain. Just work harder" (Randy Pausch).

network1.PNG

network2.PNG

network3.PNG

network0.PNG

network5.PNG

BCEM_Nadine.png

BrandonPrice1_1-1758140838545.png

CalvinHarmin_6-1758140690668.png

IMG_0494.JPG

GlitchKing_0-1606972951457.png

Anonymous User · Answer

I'm not sure I understand what you are trying to do here. I would just intersect first, then try the SymDiff.  But to do what you are trying to do you could try something like this (if you just want to print the areas):
def GetOverlap(infc_a,infc_b):
    try:
        arcpy.overwriteOutput = True
      
        # Geom object A
        Geom_a = arcpy.Geometry()
        Geom_a_list = arcpy.CopyFeatures_management(infc_a,Geom_a)

# Geom object B
        Geom_b = arcpy.Geometry()
        Geom_b_list = arcpy.CopyFeatures_management(infc_b,Geom_b)

# Find Overlaps    
        for geom in Geom_a_list:
            for geomb in Geom_b_list:
                if geom.overlaps(geomb):
                    arcpy.Intersect_analysis([geom, geomb], r'in_memory\Intersection')
                    arcpy.SymDiff_analysis(geom,geomb,r'in_memory\SymDiff')
                    Rows = arcpy.SearchCursor(r'in_memory\Intersection')
                    for Row in Rows:
                        intersectionArea = Row.Shape.area
                    Rows2 = arcpy.SearchCursor(r'in_memory\SymDiff')
                    for r in Rows2:
                        symArea = r.Shape.area
                    print "Overlap Area: %s" %str(intersectionArea/43560) # in acres
                    print "SymDiff Area: %s
"% str(symArea/43560) # in acres
        
        print 'done'
       
        
    except:
        tb = sys.exc_info()[2]
        tbinfo = traceback.format_tb(tb)[0]
        pymsg = "PYTHON ERRORS:
Traceback Info:
" + tbinfo + "
Error Info:
     " +  \
                str(sys.exc_type) + ": " + str(sys.exc_value) + "
"
        msgs = "ARCPY ERRORS:
" + arcpy.GetMessages(2) + "
"

arcpy.AddError(msgs)
        arcpy.AddError(pymsg)

print msgs
        print pymsg
        
        arcpy.AddMessage(arcpy.GetMessages(1))
        print arcpy.GetMessages(1)

if __name__ == '__main__':  # This will use the variables in this script
                
    import arcpy, sys, traceback
    arcpy.env.overwriteOutput = True

from arcpy import env
    env.workspace =r"C:\Users\Luke Kaim\Documents\University of Maine\Fall_2012\Volunteer Geographic Information\simalar"
    
    polygon = r"C:\Users\Luke Kaim\Documents\University of Maine\Fall_2012\Volunteer Geographic Information\simalar\polygon.shp"
    GNIS = r"C:\Users\Luke Kaim\Documents\University of Maine\Fall_2012\Volunteer Geographic Information\simalar\GNIS.shp"
    
    try:
        GetOverlap(polygon,GNIS)
        
    except:
        tb = sys.exc_info()[2]
        tbinfo = traceback.format_tb(tb)[0]
        pymsg = "PYTHON ERRORS:
Traceback Info:
" + tbinfo + "
Error Info:
     " +        str(sys.exc_type) + ": " + str(sys.exc_value) + "
"
        msgs = "ARCPY ERRORS:
" + arcpy.GetMessages(2) + "
"

arcpy.AddError(msgs)
        arcpy.AddError(pymsg)

print msgs
        print pymsg
        
        arcpy.AddMessage(arcpy.GetMessages(1))
        print arcpy.GetMessages(1)

Note, I displayed area in acres because the areas of my test FCs were in feet.  Also, look at the indentation.  There were a few indentation errors in yours.

Anonymous User · Answer

Original User: luke.kaim

I am still struggling with this code. Does anyone have suggestions?