2

Given a polygon with holes

from shapely.geometry import Polygon
import geopandas as gpd

# Polygon exterior:
p = [[20,767],[54,744],[107,707],
 [190,654],[265,609],[363,548],
 [462,484],[514,447],[603,389],
 [682,337],[726,310],[757,290],
 [786,277],[820,259],[843,249],
 [881,231],[921,215],[975,197],
 [1048,174],[1089,163],[1141,152],
 [1212,137],[1270,121],[1271,64],
 [1207,78],[1163,89],[1096,103],
 [1048,115],[1001,129],[949,144],
 [905,157],[874,170],[830,187],
 [781,208],[730,236],[696,255],
 [652,282],[606,306],[561,340],
 [512,370],[478,393],[436,418],
 [385,453],[330,490],[285,521],
 [229,566],[183,603],[123,652],
 [70,698],[13,749]]

# Define interior "holes":
interiors = {}
interiors[0] = [[290,543],[301,560],[393,501],[377,482]]
interiors[1] = [[507,392],[549,363],[553,367],[572,352],[588,372],[522,415]]
interiors[2] = [[599,340],[636,316],[648,334],[612,357]]
interiors[3] = [[714,262],[727,284],[821,238],[811,215]]
interiors[4] = [[850,218],[935,185],[937,187],[850,221]]
interiors[5] = [[959,159],[1066,129],[1071,146],[966,177]]
interiors[6] = [[1119,133],[1175,122],[1178,123],[1119,134]]
interiors[7] = [[1211,102],[1266,91],[1267,97],[1212,108],[1211,102]]

i_p = {k: Polygon(v) for k, v in interiors.items()}


zone = Polygon(p, [zone.exterior.coords for zone in i_p.values() \
                    if zone.within(Polygon(p)) is True])

polygon = gpd.GeoDataFrame(index=[0], crs='epsg:4326', geometry=[zone])

Extract the vertices and edges as a i) coordinate list and ii) list of indices that define the edge (end points of the line segment).

l = 0    # this will be a number ~ len(gdf)
dps = 2

all_coords = []
idx = []

for i, row in polygon.iterrows():
    exterior = row.geometry.exterior
    ring = list(exterior.coords)
    ring.pop() # pop out duplicate vertice
    if exterior.is_ccw == False:
        #-- to get proper orientation of the normals
        ring.reverse()

    for j, v in enumerate(ring):
        idx.append([l+j, l+j+1])
        rounded_x = v[0]# round(float(v[0]), dps)
        rounded_y = v[1] #round(float(v[1]), dps)
        all_coords.append([rounded_x, rounded_y])
        
    idx.append([l+j+1, l])
    l = len(all_coords) + l

    for each, inte in enumerate(row.geometry.interiors):
        ring = list(inte.coords)
        ring.pop()
        if inte.is_ccw == True:
            #-- to get proper orientation of the normals
            ring.reverse()
        
        for p in ring:
            idx.append([l+each, l+each+1])
            rounded_x = p[0] #round(float(v[0]), 2)
            rounded_y = p[1] #round(float(v[1]), 2)
            all_coords.append([rounded_x, rounded_y])#, rounded_z]) 
    
        idx.append([l+each+1, l])
        l = len(all_coords) + l
      
print(all_coords)
print('')
print(idx)

[[13.0, 749.0], [70.0, 698.0], [123.0, 652.0], [183.0, 603.0], [229.0, 566.0], [285.0, 521.0], [330.0, 490.0], [385.0, 453.0], [436.0, 418.0], [478.0, 393.0], [512.0, 370.0], [561.0, 340.0], [606.0, 306.0], [652.0, 282.0], [696.0, 255.0], [730.0, 236.0], [781.0, 208.0], [830.0, 187.0], [874.0, 170.0], [905.0, 157.0], ... [966.0, 177.0], [1071.0, 146.0], [1066.0, 129.0], [959.0, 159.0], [1119.0, 134.0], [1178.0, 123.0], [1175.0, 122.0], [1119.0, 133.0], [1212.0, 108.0], [1267.0, 97.0], [1266.0, 91.0], [1211.0, 102.0]]

[[0, 1], [1, 2], [2, 3], [3, 4], [4, 5], [5, 6], [6, 7], [7, 8], [8, 9], [9, 10], [10, 11], [11, 12], [12, 13], [13, 14], [14, 15], [15, 16], [16, 17], [17, 18], [18, 19], [19, 20], [20, 21], [21, 22], [22, 23], [23, 24], [24, 25], [25, 26], [26, 27], [27, 28], [28, 29], [29, 30], [30, 31], [31, 32], [32, 33], [33, 34], [34, 35], [35, 36], [36, 37], [37, 38], [38, 39], [39, 40], [40, 41], [41, 42], [42, 43], [43, 44], [44, 45], [45, 46], [46, 47], [47, 48], [48, 49], [49, 50], [50, 0], [50, 51], [51, 52], [52, 53], [53, 54], [54, 50], [104, 105], [105, 106], [106, 107], [107, 108], [108, 109], [109, 110], [110, 104], [164, 165], [165, 166], [166, 167], [167, 168], [168, 164], [228, 229], [229, 230], [230, 231], [231, 232], [232, 228], [296, 297], [297, 298], [298, 299], [299, 300], [300, 296], [368, 369], [369, 370], [370, 371], [371, 372], [372, 368], [444, 445], [445, 446], [446, 447], [447, 448], [448, 444], [524, 525], [525, 526], [526, 527], [527, 528], [528, 524]]

Notice print(idx) goes ....[50, 51], [51, 52], [52, 53], [53, 54], [54, 50], [104, 105], [105, 106], [106, 107], [107, 108], [108, 109], [109, 110], [110, 104]... and then again ...[167, 168], [168, 164], [228, 229], [229, 230], [230, 231], [231, 232], [232, 228], [296, 297]...

Where are the missing vertices? Why is the loop skipping?

1 Answer 1

2

You have a problem with your for loops.

As there is only one polygon:

exterior = polygon.geometry[0].exterior
interior = polygon.geometry[0].interiors
ext_coords = list(exterior.coords)[:-1] # eliminate last value

all_coords = []
idx = []
# exterior
for i,coord  in enumerate(ext_coords):
    idx.append([i,i+1])
    all_coords.append(coord)
# interior
l = len(idx)
for int in interior:
   int_coords=list(int.coords)[:-1]
   for i in int_coords:
     idx.append([l,l+1])
     l = l+1
     all_coords.append(i)
print(idx)
[[0, 1], [1, 2], [2, 3], [3, 4], [4, 5], [5, 6], [6, 7], [7, 8], [8, 9], [9, 10], [10, 11], [11, 12], [12, 13], [13, 14], [14, 15], [15, 16], [16, 17], [17, 18], [18, 19], [19, 20], [20, 21], [21, 22], [22, 23], [23, 24], [24, 25], [25, 26], [26, 27], [27, 28], [28, 29], [29, 30], [30, 31], [31, 32], [32, 33], [33, 34], [34, 35], [35, 36], [36, 37], [37, 38], [38, 39], [39, 40], [40, 41], [41, 42], [42, 43], [43, 44], [44, 45], [45, 46], [46, 47], [47, 48], [48, 49], [49, 50], [50, 51], [51, 52], [52, 53], [53, 54], [54, 55], [55, 56], [56, 57], [57, 58], [58, 59], [59, 60], [60, 61], [61, 62], [62, 63], [63, 64], [64, 65], [65, 66], [66, 67], [67, 68], [68, 69], [69, 70], [70, 71], [71, 72], [72, 73], [73, 74], [74, 75], [75, 76], [76, 77], [77, 78], [78, 79], [79, 80], [80, 81], [81, 82], [82, 83], [83, 84]]

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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