Spatial Joins for Lines & Polygons

[SaraJL]

Hello everyone!

I'm brainstorming a scenario and the best way to address it - I'm working with a professor that has a line file that needs a spatial join with a polygon layer. Easy enough, however, we're running into an issue (that I was expecting to happen) - many of the lines are overlapping with the polygon boundaries so I'm trying to reconcile that.

For example, I ran the intersect to see where the spots are where the line overlaps with the polygon boundary. A point represents every part of the line that zig zags over the polygon boundary. The polygon boundary is highlighted in pink and the line is purple:

If it was a typical situation, I would split the lines at the points, and then it would give me a little more accurate count if I ran the spatial join and where it is located with the polygon. But I know in this case, the overlapping lines are a bit of a sticking point.

Any recommendations? I need to create some type of join between the lines and polygons, but there isn't a common field in either dataset. There are also many lines over one polygon, and one line can stretch multiple polygons. Ideally, I need some kind of table where it's one polygon (census tract in this case) with an entry for each line. Or just some why to create a tie between the layers when they only share a spatial relationship and no shared attributes.

To be honest - I'm not sure if there is much I can do about it, there is probably going to be some type of margin of error, but just exploring possibilities.

Also, something that does not involved writing a script please!

Thank you!

[D_Atkins]

There's nothing wrong with spatial joins between lines and polygons, other than some additional decision making.

We perform similar analyses routinely (think of Residential Waste Service Routes overlaid on a street network); so I had the work ready to go -- the best way to share is a minimum-working-example by script, copy-and-paste, turnkey ready for a Notebook:

The first 40 lines are simply preprocessing: I take census block groups and secondary roads from the Living Atlas; rather than process the entire dataset, I make view-layers, apply a selection, then copy the selected features to a local workspace.  You can replace all of this by referencing your own features in lines 46/47.

Lines 46 to 60 'do the work' with just two tools: (Intersect, then, Spatial Join).  I've attached some screenshots to follow, which illustrate the issue with your client's request: the highway follows the polygon boundaries, gets split almost at random, then gets attributes joined from 'one' polygon when really it roughly belongs to both.

import arcpy

poly = "https://services2.arcgis.com/FiaPA4ga0iQKduv3/arcgis/rest/services/USA_Block_Groups_v1/FeatureServer/0"
line = "https://services2.arcgis.com/FiaPA4ga0iQKduv3/arcgis/rest/services/Transportation_v1/FeatureServer/3"


### PREPROCCING FOR LIVING ATLAS (!

  ## make a layer from census blocks:
poly_layer = arcpy.management.MakeFeatureLayer(poly, "Atlas Census Blocks LayerView")

  ## select from that layer:
poly_selection = arcpy.management.SelectLayerByAttribute(
    in_layer_or_view = poly_layer, 
    selection_type   = "NEW_SELECTION",
    where_clause     = "COUNTY = '173' And STATE = '37'"
)

  ## copy selection to memory:

  # select returns a 'updated layer' as var[0], and 'count' as var[1]
print("Selected:", poly_selection[1], "counties")

  # use the selected layer to select all roads in the same area:
line_layer = arcpy.management.MakeFeatureLayer(line, "Atlas Roads LayerView")

  ## select from that layer:
line_selection = arcpy.management.SelectLayerByLocation(
    in_layer         = line_layer, 
    overlap_type     = "INTERSECT",
    select_features  = line_layer,
    selection_type   = "NEW_SELECTION"
)

print("Selected:", line_selection[2], "roads")

arcpy.management.CopyFeatures(poly_selection[0], "memory/my_poly")
arcpy.management.CopyFeatures(line_selection[0], "memory/my_lines")

### FINISHED PREPROCESSING: TWO LAYERS WITH SELECTIONS FOR INPUT...


### SPLIT LINES BY POLYGONS:

## References to lines and polygons:
my_lines = "memory/my_lines"
my_poly  = "memory/my_poly"


  ## split line_layer by poly_layer, create new 'lines_split' featureclass.
split_lines = arcpy.analysis.Intersect([my_lines, my_poly], "memory/lines_split")

  ## join the polygon attributes to the lines;
join_attr = arcpy.analysis.SpatialJoin(
    target_features = split_lines, 
    join_features   = poly_layer,
    out_feature_class = "memory/lines_split_wPolyAttr",
    join_operation = "JOIN_ONE_TO_ONE",
    join_type      = "KEEP_ALL",
    match_option   = "INTERSECT",
    search_radius  = "0 Feet"
)
    
print("Completed")    
      

Symbolized results and full table of contents:

 

The tools appear to do what you're asking for.  The results are less than ideal, of course.

I would perhaps suggest buffering road lines and using these road-buffer polygons for further analysis; this is more reflective of the real world of roads (which generally include 50-100 feet of right-of-way, and would settle many of the questions raised ... at the expense of geoprocessing complexity).