Creating polygon from 2 txt files (Parent ID + Coordinates) using ArcPy?

Question

I want to create polygons based on two separate text files:

• the first one (as parents) lists the Polygon ID (ID_el) and three to four Points ID (IDpt1, IDpt2, IDpt3..)
• the second contains the Points ID and their xyz coordinates (RW, HW,HH).

I want for each ID in the first file to set an Array (3 to 4 rows) of the coordinates in the second file

I found different scripts to do some of the parts of the process, but I really need help to figure out how to get the Points ID and coordinates for each polygon into an Array, keeping also the attributes of the parent element (ElTyp, ID_el).

Here is the bit of code I used to extract these 2 text files:

with open(arcpy.GetParameterAsText (0),'r') as f1:
    with open(r'Nodes.txt', 'a') as f2:
        f2.write("ID RW HW HH" + "\n")
        for line in f1:
            if line.startswith("ND"):
                f2.write(line.replace("ND", "").lstrip())
with open (arcpy.GetParameterAsText (0),'r') as f1:
    with open (r'Polygonen.txt', 'a') as f3:
        f3.write("ELtyp ID_el IDpt1 IDpt2 IDpt3 IDpt4 Mat" + "\n")
        for line in f1:
            coord = " ".join(line.split() [0:-1])
            mat = " ".join(line.split() [-1:])
            if line.startswith ("E3T"):
                f3.write ((coord.lstrip())+ "  " + mat+ "\n")
            elif line.startswith ("E4Q"):
                f3.write(line.lstrip())

Answer 1

A similar approach to what has Erica published, but with more details:

1. Export text files to a file geodatabase table (useful because you will get 0 for the vertex ID for those polygons which have only 3 points - thus no need to handle this later on);
2. Convert those tables into dictionaries;

polygons {41880: (26287, 26286, 26748, 26747), 41879: (26748, 26286, 26747, 0)}

vertices {26747: (19, 20, 21), 26748: (12, 14, 15), 26286: (13, 16, 17), 26287: (22, 24, 25)}

3. Building an array of polygons (each has points with XYZ values);
4. Pre-create the polygon feature class with the Z coordinates stored (important); doable with the GP tool and adding the PolyID field.
5. When loading polygons, important to specify the has_z attribute, otherwise it won't store Z-coordinates. You can see the Z-value for every vertex of the polygon when in the Editing session and having the Edit Sketch Properties window open (while editing the polygon with vertices shown);
6. As the last step, use the Join Field GP tool to transfer all other attribute fields from the Polygons file GDB table (based on the PolyID). I didn't want to take them into the dictionary just to keep things more clear.

Useful Esri Help links: Polygon class; Writing geometries; da.Insert cursor

The ready-to-use code:

import arcpy
import os
folder = "C:\GIS\Temp"
os.chdir(folder)
polygonsFile = "polygons.txt"
verticesFile = "vertices.txt"

outfileGDB = r"C:\GIS\Temp\test.gdb"

if not arcpy.Exists(r"C:\GIS\Temp\test.gdb\polygons"):
    arcpy.TableToTable_conversion(in_rows=polygonsFile,
                                  out_path="C:\GIS\Temp\test.gdb",
                                  out_name="polygons",
                                  where_clause="")

    arcpy.TableToTable_conversion(in_rows=verticesFile,
                                  out_path="C:\GIS\Temp\test.gdb",
                                  out_name="vertices",
                                  where_clause="")

arcpy.env.workspace = outfileGDB
with arcpy.da.SearchCursor("polygons","*") as poly_cur:
    print "polygons"
    poly_dict = {x[1]: x[2:] for x in poly_cur}
    print poly_dict

with arcpy.da.SearchCursor("vertices","*") as vertex_cur:
    print "vertices"
    vertex_dict = {x[1]: x[2:] for x in vertex_cur}
    print vertex_dict

polyArray = {}
for polykey in poly_dict:
    polynodes = poly_dict[polykey]
    coordsList = []
    for polynode in polynodes:
        coords = (v for k,v in vertex_dict.iteritems() if k == polynode)
        for coord in coords:
            coordsList.append(coord)
    polyArray[polykey] = coordsList

print polyArray
#>>>polyArray
#>>>{41880: [(22, 24, 25), (13, 16, 17),(12, 14, 15), (19, 20, 21)],
#>>>41879: [(12, 14, 15), (13, 16, 17), (19, 20, 21)]}

features = []
fc = r"C:\GIS\Temp\test.gdb\PolygonFC"

for key,values in polyArray.iteritems():
    features.append((key,arcpy.Polygon(arcpy.Array([arcpy.Point(*value) for value in values]),
                                       arcpy.SpatialReference(4326),True))) #True - for has_z parameter

with arcpy.da.InsertCursor(fc,["PolyID","SHAPE@"]) as cur:
    for feature in features:
        cur.insertRow(feature)

Answer 2

This is an outline/pseudocode of one possible method...

1. Create dictionaries for the point ID and coordinates. Populate as you loop through the table.

   dictRW = {}
   dictHW = {}
   dictHH = {}
   
   for row in arcpy.da.SearchCursor(tblCoordinates):
       ID = row[0]
       dictRW[ID] = row[1]
       dictHW[ID] = row[2]
       dictHH[ID] = row[3]
   

2. Create dictionaries using the polygon ID and its associated point IDs. Populate as you loop through the table.

   listPolygons = []
   ELType = {}
   pt1 = {}
   pt2 = {}
   pt3 = {}
   pt4 = {}
   
   for row in arcpy.da.SearchCursor(tblPolygon):
       polyID = row[1]
       ELType[polyID] = row[0]
       pt1[polyID] = row[2]
       pt2[polyID] = row[3]
       pt3[polyID] = row[4]
       pt4[polyID] = row[5]
       listPolygons.append(polyID)
   

3. Work through the list of polygons and points to insert the data into a feature class. (The following assumes the feature class exists and has the fields that you want. If you want to create it in-script, refer to arcpy.CreateFeatureclass_management)

   for polyID in listPolygons:
       featureData = []
   
       featureData.append(polyID) # first attribute is the polygon ID
   
       ptID1 = pt1[polyID]
       ptID2 = pt2[polyID]
       ptID3 = pt3[polyID]
       ptID4 = pt4[polyID]
   
       if not ptID4 is None: # this is an attempt to deal with 3-point polygons
           array = arcpy.Array([arcpy.Point(dictRW[ptID1], dictHW[ptID1], dictHH[ptID1]),
                                arcpy.Point(dictRW[ptID2], dictHW[ptID2], dictHH[ptID2]),
                                arcpy.Point(dictRW[ptID3], dictHW[ptID3], dictHH[ptID3]),
                                arcpy.Point(dictRW[ptID4], dictHW[ptID4], dictHH[ptID4])])
       else:
           array = arcpy.Array([arcpy.Point(dictRW[ptID1], dictHW[ptID1], dictHH[ptID1]),
                                arcpy.Point(dictRW[ptID2], dictHW[ptID2], dictHH[ptID2]),
                                arcpy.Point(dictRW[ptID3], dictHW[ptID3], dictHH[ptID3])])
   
       polygon = arcpy.Polygon(array)
       featureData.append(polygon)
   
       with arcpy.da.InsertCursor(outFC, ["PolyID", "SHAPE@"]) as cursor:
           cursor.insertRow(featureData)
   

Note that the InsertCursor fields should be adapted if you want to preserve/insert additional attributes (e.g., ELType). Just make sure you also include a featureData.append() statement in order -- that is building the list of data which the insertRow command uses.

Answer 3

As my understanding of your problem, I wrote a script below. I assumed that you have X,Y and Z values for every point.

Script:

# Example of PrentID text file
#ELtyp ID_el IDpt1 IDpt2 IDpt3
#E3T 41879 26748 26286 26747
#E3T 41880 26287 26286 26748
#E3T 41881 26749 26287 26748
#E3T 41882 26749 26750 26287

#-----------------------------------------------------------
# Example of PointsID text file
#ID RW HW HH
#26748 4.40384415e+006 5.79104121e+006 7.85830000e+001
#26749 4.40386009e+006 5.79106456e+006 7.83740000e+001
#26750 4.40386214e+006 5.79108649e+006 7.80960000e+001
#26286 4.40381493e+006 5.79104349e+006 7.85940000e+001
#26287 4.40383655e+006 5.79107131e+006 7.84350000e+001
#26747 4.40383190e+006 5.79101331e+006 7.87290000e+001

#-----------------------------------------------------------
import arcpy

arcpy.env.overwriteOutput = True
def main():
    # Polygon ID text file path
    parent = r"D:\Python\Scratch\ParentID.txt"
    # Points ID text file path
    point = r"D:\Python\Scratch\Points.txt"
    # spatial reference wkid
    wkid = 31468
    # output path (could be gdb or folder which creates
    # feature class or shapefile respectively)
    output = r"D:\Python\ScratchDatabase\Geodatabase.gdb\polygon_123"


    points = {}
    fields = []
    # create point id and coordinates dictionary
    with open(point, "rt") as point_file:
        for line in point_file:
            if line.upper().startswith('ID'.upper()):
                continue
            points.update({line.split()[0]: [float(pt) for pt in line.split()[1:]]})
    sr = arcpy.SpatialReference(wkid)
    ftr = arcpy.CreateFeatureclass_management('in_memory', 'Test', "POLYGON", spatial_reference=sr)
    flds = [fld.name for fld in arcpy.ListFields(ftr)]
    with open(parent, 'rt') as parent_file:
        for line in parent_file:
            if line.upper().startswith("ELtyp".upper()):
                for field in line.split():
                    if field not in flds:
                        # Add Fields for attributes
                        arcpy.AddField_management(ftr, field, "TEXT", field_length=10, field_alias=field)
                    fields.append(field)
                fields.append("SHAPE@")
            else:
                # insert values
                with arcpy.da.InsertCursor(ftr, fields) as cursor:
                    row = [val for val in line.split()]
                    pts = [arcpy.Point(points.get(val)[0],
                                       points.get(val)[1],
                                       points.get(val)[2]) for val in line.split()[2:]]
                    row.append(arcpy.Polygon(arcpy.Array(pts), sr, True))
                    cursor.insertRow(row)
    arcpy.CopyFeatures_management(ftr, output)
    print "Done"


if __name__ == '__main__':
    main()

May be this will help.