4

I have a shapefile with polygons representing unique populations. I want to determine contiguous polygons in my shapefile, and then extract the border separating contiguous polygons as a separate polyline file. Importantly, I require the created polyline file to indicate the two polygons it separates.

I am aware of Polygon Neighbor to produce a table of all contiguous pairs, but this doesn't help me extract borders as polylines.

Consider group A in the following example map:

Example Border

Group A has 8 contiguous pairs: AB, AC, AD, AE, AF, AG, AH, AI. I would like the output file to separate the aggregate border for group A into 8 separate lines, each line referencing only the shared border between each contiguous pair; e.g., the border of AH in red would be a single line in the output data. As noted before it is important that the output file also identifies the two groups each line separates.

To clarify, I require a general solution for all possible contiguous pairs in my dataset. I do not have a single group such as the example group A that I am interested in.

This question is similar to Extracting polyline of border between differing polygons using R or GRASS?, except I'm interested in a solution using ArcGIS 10.2 for Desktop. I have some knowledge of Python if an ArcPy solution is easiest.

Community
  • 1
acd
  • 107
  • 6
  • 1
    If you need to identify from the polylines which polygons are on either side use Polygon to Line http://help.arcgis.com/EN/arcgisdesktop/10.0/help/index.html#//00170000003t000000 and specify IDENTIFY_NEIGHBORS (advanced license required) It would be best if you didn't use shapefiles as the FIDs are not static, preferably convert the polygon shapefile to a geodatabase feature class (with static OIDs) so that the ids are reliable. – Michael Stimson Sep 08 '15 at 22:12
  • @MichaelMiles-Stimson: If I am able to figure out how to extract borders as polylines I will follow your suggestion. Thanks! – acd Sep 08 '15 at 22:39
  • That's what the tool does. Where polygons coalesce only one line is produced and if you specify it will store the OID/FID of the source polygons in a field. – Michael Stimson Sep 08 '15 at 22:48
  • @MichaelMiles-Stimson: Sorry I didn't read your comment correctly the first time around. Thank you for the response, I was able to almost do what I want with your suggestion. I say almost because approximately 1/5 of the polylines produced are disjoint; e.g., a single border segment might be split into 3 disjoint polylines, where one line has the appropriate left/right FIDs, but then the other two have only one of the FIDs and an unmatched FID of -1 for the other. Is this a common problem you are aware of? – acd Sep 09 '15 at 20:13
  • Yes it's a common problem, because of tiny gaps in the polygons the lines diverge and are populated with -1 (no polygon) on the outside. A tool like Eliminate or Integrate (see http://gis.stackexchange.com/questions/11004/removing-small-spaces-slivers-between-polygons or http://gis.stackexchange.com/questions/154264/how-to-correct-topology-area-boundary-errors-slivers-gaps-on-arcs-in-arcgis/154309#154309 - read the warnings) will help close the gaps, after the gaps are removed you should have single lines at the boundaries. – Michael Stimson Sep 09 '15 at 21:39
  • @MichaelMiles-Stimson: thank you for the assistance, I will take a look at Eliminate and Integrate. – acd Sep 10 '15 at 14:04

1 Answers1

3

This will convert your polygons into arcs very quickly. You need 2 layers in table of content: 1st - polylines, 2nd - polygons

import arcpy, traceback, os, sys,time
from arcpy import env
env.overwriteOutput = True
dissolved="in_memory/dissolved"

try:
    def showPyMessage():
        arcpy.AddMessage(str(time.ctime()) + " - " + message)

    mxd = arcpy.mapping.MapDocument("CURRENT")
    theLinksLayer=arcpy.mapping.ListLayers(mxd)[0]
    PGONS=arcpy.mapping.ListLayers(mxd)[1]
    arcpy.DeleteFeatures_management(theLinksLayer)
    curT = arcpy.da.InsertCursor(theLinksLayer,"SHAPE@")
    with arcpy.da.SearchCursor(PGONS, "Shape@") as cursor:
        for row in cursor:
            feat = row[0].boundary()
            curT.insertRow((feat,))
    arcpy.Dissolve_management(theLinksLayer, dissolved,"","", "SINGLE_PART")
    arcpy.DeleteFeatures_management(theLinksLayer)
    m=0
    arcpy.AddMessage("Dissolving ...")
    curT = arcpy.da.InsertCursor(theLinksLayer,"SHAPE@")
    with arcpy.da.SearchCursor(dissolved, "SHAPE@") as cursor:
        for row in cursor:
            m+=1
            curT.insertRow((row[0],))

    arcpy.Delete_management("in_memory")
    arcpy.RefreshActiveView()
    arcpy.AddMessage("\nPolygons unbuilt into %i lines\n" %m)
    del curT

except:
    message = "\n*** PYTHON ERRORS *** "; showPyMessage()
    message = "Python Traceback Info: " + traceback.format_tb(sys.exc_info()[2])[0]; showPyMessage()
    message = "Python Error Info: " +  str(sys.exc_type)+ ": " + str(sys.exc_value) + "\n"; showPyMessage()

From here it is enough to spatially join lines to polygons (one to many, SHARE_A_LINE_SEGMENT_WITH) and dissolve output using line name and first, last polygon IDs to get this:

enter image description here

FelixIP
  • 22,922
  • 3
  • 29
  • 61
  • Thank you for the very thorough response. I apologize but I am not following something very basic in your explanation. I am confused where the 1st polyline layer in the Table of Contents is coming from? My aim is to create a polyline layer of contiguous border segments so I don’t understand how to add it to the Table of Contents in the first place. Thank you. – acd Sep 09 '15 at 00:40
  • Create empty polyline feature layer. Make it sitting on the top of all the other layers. Place your source polygons layer below it. – FelixIP Sep 09 '15 at 01:06
  • To be honest there many other ways to derive polygon edges as Michael said. You might start editing empty polyline layer. Select all polygons, copy-paste to polyline layer. Go advanced editing, pick planarise. I offered script simply because it is sitting in one of my custom toolboxes and I preffer to use it – FelixIP Sep 09 '15 at 01:18
  • Thank you. I am learning to use script so I am slow but I prefer to do it this way. I created an empty polyline feature layer, placed it at the top, followed by my source polygon layer. When running the script I keep getting the error message: ERROR 999999: Error executing function. Cannot acquire a lock. Cannot acquire a lock. Failed to execute (Dissolve). I am unable to determine why the layers are locked when using the Dissolve tool? – acd Sep 09 '15 at 01:56
  • Try to change dissolved="in_memory/dissolved" to something like dissolved=r"c:\work\dissolved.shp". in_memory storage can be tricky. Also disable arcpy.Delete_management(in_memory). Note that I've added escept block at the end of script. No idea why it wasn't best at the start. Python is good, don't give up – FelixIP Sep 09 '15 at 02:00
  • Except block noted thanks. I've spent a good amount of time trying your suggestion with various other attempts to troubleshoot but I still get the same error message regardless. I noticed that prior to the previously stated error, I'm getting the following from the DeleteFeatures tool: WARNING 000117: Warning empty output generated. The result following this error is that the empty polyline feature layer now maps the outline of my shapefile polygons, and the attribute table has exactly the same number of rows as the original shapefile and includes a variable Id equal to 0 for all rows. – acd Sep 09 '15 at 03:11
  • Check environment extent in geoprocessing. Set it to union of inputs – FelixIP Sep 09 '15 at 03:13
  • I really wish I had better news than the extent was set to default so I changed it to union of inputs, but I'm still getting the same warning and error message as before. This is true whether I use your original code, or your suggestion to change dissolved="in_memory/dissolved" and disable arcpy.Delete_management(in_memory). I even tried closing and reopening ArcMap to be sure no layers were locked before executing the script. – acd Sep 09 '15 at 03:43
  • Let us continue this discussion in chat. – FelixIP Sep 09 '15 at 03:44