4

I wondered if it is possible to determine the number of sides / edges (not vertices) that a polygon has. Imagine the following case. There is a shapefile with 10 polygons.

Is there any way to determine for each polygon, how many sides there are and then write these into the attribute table?

I would be happy for any solution using open source tools.

In addition, I aim to use that in for describing voronoi polyons.

  • 8
    if the polygon has no hole, it is simply number of vertices -1 – gene Oct 12 '14 at 19:16
  • 5
    The general form of that formula is nsides = nverts - nrings (works for both single and multi-part polygons, with and without holes) – Vince Oct 13 '14 at 0:06
5

In QGIS if your polygons do not have any holes or multi-parts:

l = iface.activeLayer()
for f in l.getFeatures():
    print f['NAME']
    print 'no. edges: %d' %(len(f.geometry().asPolygon()[0])-1)

replace 'NAME' with some identifier in your layer attribute table.

Concerning writing to the attribute table check the instructions in the PyQGIS Cookbook - Modifying Vector Layers.

6

With Python and Fiona, Polygons and MutiPolygons (multi-parts) are different geometries:

1) multi-parts geometries

 import fiona
 shape = fiona.open("polygons.shp")
 # shapefile schema
 print c.schema
 {'geometry': 'Polygon', 'properties': OrderedDict([(u'id', 'int:10')])}
 # first feature
 first = shape.next()
 print first (GeoJSON format)
 {'geometry': {'type': 'Polygon', 'coordinates': [[(244697.45179524383, 1000369.2307574936), (244827.15493968062, 1000373.0455558595), (244933.96929392271, 1000353.9715640305), (244933.96929392271, 1000353.9715640305), (244930.15449555693, 1000147.9724522779), (244697.45179524383, 1000159.4168473752), (244697.45179524383, 1000369.2307574936)]]}, 'type': 'Feature', 'id': '0', 'properties': OrderedDict([(u'id', 1)])}


 multi = fiona.open("multipolygons.shp")
 print c.schema
 {'geometry': 'Polygon', 'properties': OrderedDict([(u'id', 'int:10')])}
 # first feature
 first = shape.next()
 print first
 {'geometry': {'type': 'MultiPolygon', 'coordinates': [[[(244697.45179524383, 1000369.2307574936), (244827.15493968062, 1000373.0455558595), (244933.96929392271, 1000353.9715640305), (244933.96929392271, 1000353.9715640305), (244930.15449555693, 1000147.9724522779), (244697.45179524383, 1000159.4168473752), (244697.45179524383, 1000369.2307574936)]], [[(246082.22360202507, 1000453.1563215409), (246139.44557751188, 1000460.7859182726), (246189.03795626713, 1000403.5639427857), (246189.03795626713, 1000403.5639427857), (246086.03840039085, 1000132.7132588148), (245990.66844124615, 1000205.1944277647), (246082.22360202507, 1000453.1563215409)]]]}, 'type': 'Feature', 'id': '0', 'properties': OrderedDict([(u'id', 1)])}

2) With simple Polygons, the coordinates are

coord = first['geometry']['coordinates']
[[(244697.45179524383, 1000369.2307574936), (244827.15493968062, 1000373.0455558595), (244933.96929392271, 1000353.9715640305), (244933.96929392271, 1000353.9715640305), (244930.15449555693, 1000147.9724522779), (244697.45179524383, 1000159.4168473752), (244697.45179524383, 1000369.2307574936)]
print len(coord)
1 # -> one polygon

And the equivalent LinearRing coordinates:

linearR = coor[0] # or coord = first['geometry']['coordinates'][0]
print linearR
[(244697.45179524383, 1000369.2307574936), (244827.15493968062, 1000373.0455558595), (244933.96929392271, 1000353.9715640305), (244933.96929392271, 1000353.9715640305), (244930.15449555693, 1000147.9724522779), (244697.45179524383, 1000159.4168473752), (244697.45179524383, 1000369.2307574936)]

So

nb_vertices = len(lineaR) #(x,y) points
print nb_vertices
7 
nb_edges = nb_vertices -1
print nb_edges
6

3) With MultiPolygons, use a for loop

4) If the Polygons have holes, use Shapely with Fiona

holes = fiona.open("poly_holes.shp")
from shapely.geometry import shape
# First feature
first = holes.next()
print first['geometry']['coordinates']
# conversion to shapely geometry
shape = first['geometry']
[[(1.0, 1.0), (1.0, 7.0), (7.0, 7.0), (7.0, 1.0), (1.0, 1.0)], [(2.0, 3.0), (4.0, 3.0), (4.0, 5.0), (2.0, 5.0), (2.0, 3.0)], [(5.0, 5.0), (6.0, 5.0), (6.0, 6.0), (5.0, 6.0), (5.0, 5.0)]]
# exterior coordinates
print list(shape.exterior.coords) # = LinearRing
[(1.0, 1.0), (1.0, 7.0), (7.0, 7.0), (7.0, 1.0), (1.0, 1.0)]
edges_ext= len(list(shape.exterior.coords))
print edges_ext
5

Same with interior coordinates (look at the Shapely manual)

  • Thank you for pointing at Fiona and Shapely, I did not know either of them. Will check it out later today. – Jens Oct 13 '14 at 10:42
  • Okay, I finally got Fiona installed and running. Thank god. Then I reproduced your example code and got a python script running that showed - number of edges for the first polygon. I guess I have to loop through all other polygons, yet I have no idea how. Hence, the accepted answer goes to underdark - tadaaa. – Jens Oct 16 '14 at 10:24
  • 1
    Why tadaaa ? It is not a competition. If you are more comfortable with PyQGIS, use the underdark solution. If you want to use only Python, you need to learn first how to use Fiona and shapely. – gene Oct 16 '14 at 13:10
  • Sorry for the tadaaa. I was actually pretty happy to get the python setup finally running and will definitely continue with python only solutions. Thx again – Jens Oct 16 '14 at 14:18
1

I finally found also an R-solution for my problem. The rgeos package offers the get.pts function which allows to extract the number of vertices. As we have learned from the first comments, the number of edges is equal to the number of vertices - 1. Programming a solution for my task in R seems now much easier.

  • Rgeos use the GEOS library, as shapely with len(list(shape.exterior.coords)(= get.pts), PyQGIS, GeoDjango, PostGIS and others... – gene Oct 16 '14 at 15:10
  • Yep, that's how I found rgeos. – Jens Oct 16 '14 at 17:47

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.