-1

I am trying to speed up some current code that I have for looping through lines, getting their segments and ultimately the vertex points in Python.

Can anyone recommend a faster method than what I currently have (perhaps using a map function or dictionary)?

My code is very slow, especially when dealing with large numbers of features.

    sup_list = [row for row in arcpy.da.SearchCursor(sup_lyr, fields)]

    # with arcpy.da.SearchCursor(sup_lyr, fields) as sup_cursor:
    for sup_row in pbar(sup_list, total=result, ascii=True):
        exp = str(tracking_field) + " = {0}".format(sup_row[3])
        sup_sel = arcpy.SelectLayerByAttribute_management(sup_lyr, "NEW_SELECTION", exp)
        inf_sel = arcpy.SelectLayerByLocation_management(inf_lyr, "WITHIN_A_DISTANCE_GEODESIC", sup_sel, "50 Meters"
                                                     , "NEW_SELECTION", "NOT_INVERT")
        sup_vertexList = []
        sup_proj_geom = sup_row[1].projectAs(calc_utm(sup_row[1].firstPoint))
        for part in sup_proj_geom:
            for pnt in part:
                sup_vertexList.append([pnt.Y, pnt.X])
        with arcpy.da.SearchCursor(inf_sel, fields) as inf_cursor:
            for inf_row in inf_cursor:
                inf_vertexList = []
                inf_proj_geom = inf_row[1].projectAs(calc_utm(inf_row[1].firstPoint))
                for part in inf_proj_geom:
                    for pnt in part:
                        inf_vertexList.append([pnt.Y, pnt.X])
  • Did you try the method explained here? gis.stackexchange.com/questions/26369/… – Richard Morgan Jun 12 at 17:19
  • 1
    Show the lines prior to the code you added, are you using a cursor inside a cursor? What does the function calc_utm do? – BERA Jun 12 at 17:19
  • I have updated the code in my sample to add more at the beginning, if that is helpful. – AttilaTheHun Jun 12 at 17:35
  • 1
    A Spatial Join using the Many To One option would probably make the whole thing faster with the layer of the inner loop being the target layer and the layer with the attribute selection criteria being the join layer. Then there would be a common attribute selection on the Spatial Join features that would match your outer loop layer and you could use a dictionary to match them. The spatial join and attribute dictionary matching would run hundreds of times faster than large sets of spatial selections, and it becomes increasingly efficient as you add more features. – Richard Fairhurst Jun 12 at 17:49
  • 2
    Using Desktop, the easiest way to get a code snippet for a tool is to open the tool like you normally would (Analysis Tools.tbx->Overlay-Spatial Join), set up the parameters the way you want (Target layer, Join Layer, Output Name, JOIN_MANY_TO_ONE, optionally keep only matched features for performance, adjust field list to just include fields you need for performance, WITHIN_A_DISTANCE_GEODESIC, 50 Meters) and run it to test. Then you can open the Results tab (Geoprocessing menu->Results), right click the tool and choose Copy As Python Snippet. Paste into Python and customize as needed. – Richard Fairhurst Jun 12 at 18:34
0

Here is a version of a function I regularly use:

def get_point_coords_from_lines(shp_in, target_epsg, line_id_fieldname, output_location):
    """
    Get vertices for lines with NODE_ID order field
    :shp_in:            input polyline shapefile
    :target_epsg:       desired EPSG for point coordinates
    :line_id_fieldname: fieldname for line identifier; must be unique, like a primary key
    :output_location:   folder path for output files
    """

    import arcpy

    # copy input SHP just in case...
    temp1 = output_location + '\\' + 'temp_lines_in.shp'
    arcpy.CopyFeatures_management(shp_in, temp1, config_keyword = "")

    temp2 = output_location + '\\' + 'temp_lines_in_proj.shp'
    arcpy.Project_management(temp1, temp2, arcpy.SpatialReference(target_epsg))
    arcpy.Delete_management(temp1)

    my_pts = output_location + '\\' + 'temp_proj_as_Points.shp'
    arcpy.FeatureVerticesToPoints_management(temp2, my_pts, point_location = "ALL")
    arcpy.Delete_management(temp2)
    arcpy.AddXY_management(my_pts)

    ### this gives a "node_id" aka a point order so that you can rebuild the 
    ### lines with this as the sort field and 'line_id_fieldname' as the group field
    arcpy.AddField_management(my_pts, "NODE_ID", "SHORT", field_is_nullable="NULLABLE")
    line_id_previous = arcpy.da.SearchCursor(my_pts, (line_id_fieldname,)).next()[0]
    counter = 0
    with arcpy.da.UpdateCursor(my_pts, ['OID@', line_id_fieldname, 'NODE_ID']) as cursor:
        for row in cursor:
            line_id_current = row[1]
            if line_id_current == line_id_previous:
                counter += 1
            else:
                counter = 1
            row[2] = counter
            line_id_previous = line_id_current
            cursor.updateRow(row)

    return my_pts

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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