How to automate performing analysis (intersection) between objects within different feature classes, for large number of objects?

Question

I have two tables of roads - the same roads, from different sources, and am trying to compare how similar their positions are. I'm wanting to buffer each road in one of the feature classes and calculate the % overlap between it and the corresponding road from the other feature class. Could obviously do this manually but since there are hundreds of roads would rather use Python. The thing is, I'm new to it. Any help on how to automate using the intersect tool between objects from different tables would be very welcome. The objects will have an attribute in common.

Answer 1

I am assuming you are using ArcGIS 10.1 with this answer. I also assume that there is a 1:1 relationship between the two feature classes, that they are both in a geodatabase, that they have the same spatial reference, and that they are already sorted in the same order. I assume that if the buffers of two pairs of features do not intersect at all that you do not care about that information being output. If any of this is not true then some adjustments will need to be made.

This script will create a pair of buffered polygon feature classes as well as an intersection feature class and will attribute it with the percent overlap between pairs of features as well as FIDs for joining back to the buffered polygon feature classes (or your original lines if they share the same FIDs). The pairing of features is only made by simultaneously iterating over both of the buffered polygon feature classes (using itertools.izip) and assuming that they are pre-sorted and 1:1. I do have an assertion that will result in an error message if the values of the common key field you mentioned are not equal.

The if __name__ == "__main__": section is where you configure the script using your own data. This is a common way to test functions in a Python module.

import arcpy, os, itertools

def computeIntersectOverlapPairs(feature_class_1, key_field_1, feature_class_2, key_field_2, output_feature_class, buffer_distance):
    desc = arcpy.Describe(feature_class_1)
    buffer1 = arcpy.Buffer_analysis(feature_class_1, os.path.join(os.path.dirname(output_feature_class), "{0}_Buffer".format(os.path.basename(feature_class_1))), buffer_distance)
    buffer2 = arcpy.Buffer_analysis(feature_class_2, os.path.join(os.path.dirname(output_feature_class), "{0}_Buffer".format(os.path.basename(feature_class_2))), buffer_distance)
    createOutputFC(output_feature_class, desc.SpatialReference)
    with arcpy.da.SearchCursor(buffer1, ["OID@", "SHAPE@", key_field_1]) as rows1:
        with arcpy.da.SearchCursor(buffer2, ["OID@", "SHAPE@", key_field_2]) as rows2:
            with arcpy.da.InsertCursor(output_feature_class, ["SHAPE@", "FID_1", "FID_2", "OVERLAP"]) as outrows:
                for row1, row2 in itertools.izip(rows1, rows2):
                    assert row1[2] == row2[2], "Key values must be equal! (FID1: {0}, FID2 {1}, Key1: {2}, Key2: {3})".format(row1[0], row2[0], row1[2], row2[2])
                    shape1 = row1[1]
                    shape2 = row2[1]
                    if not shape1.disjoint(shape2):
                        intersect_shape = shape1.intersect(shape2, 4)
                        percent_overlap = (intersect_shape.area / (shape1.area + shape2.area)) * 100
                        outrows.insertRow((intersect_shape, row1[0], row2[0], percent_overlap))

def createOutputFC(path, spatial_reference=None):
    arcpy.CreateFeatureclass_management(os.path.dirname(path), os.path.basename(path), "POLYGON", spatial_reference=spatial_reference)
    arcpy.AddField_management(path, "FID_1", "LONG")
    arcpy.AddField_management(path, "FID_2", "LONG")
    arcpy.AddField_management(path, "OVERLAP", "DOUBLE")

if __name__ == '__main__':
    arcpy.env.overwriteOutput = True
    inputfc1 = r"C:\GISData\test8.gdb\atlantic_hurricanes_2000_lines"
    inputfc2 = r"C:\GISData\test8.gdb\atlantic_hurricanes_2000_lines_shifted"
    outputfc = r"C:\GISData\test8.gdb\atlantic_hurricanes_2000_lines_overlap"
    field1 = "EVENTID"
    field2 = "EVENTID"
    bufferdist = 0.1
    computeIntersectOverlapPairs(inputfc1, field1, inputfc2, field2, outputfc, bufferdist)

Visual analysis

Input Lines:

Buffered Lines:

Overlapping Area (~12% in this case):

Answer 2

You could:

1. calculate buffers for the two shapes (choose "dissolve field" as NONE, so every feature/road gets its own buffer...
2. intersect the two buffer-shapes
3. add field in resulting shapes attribute table
4. calculate the area of this shape
5. use the calculated area to create an overlap-percentage ranking (the smaller the area, the higher the overlap)