Count the number of vertices in each polygon in a shape file using Python and OGR

Someone asked me at work how to find the number of vertices in each polygon in series of shape files. I spent a few minutes to write a python script to do the job.

At first I just looped over each feature, got the geometry reference and then tried to use GetPointCount() to give me the number of points. Imagine my surprise when it returned zero. Here is my first try:

>>> from osgeo import ogr
>>> driver = ogr.GetDriverByName('ESRI Shapefile')
>>> hinputfile = driver.Open(r'c:\GeoData\AZ.shp', 0)
>>> in_layer = hinputfile.GetLayer(0)
>>> in_layer.GetFeatureCount()
5
>>> feature = in_layer.GetFeature(0)
>>> geometry = feature.GetGeometryRef()
>>> geometry.GetPointCount()
0

It turned out to be pretty simple - I had to loop over each linear ring that the polygon contains and count its points. Abbreviated, here is what I ended up with:

from osgeo import ogr

def process_file(input_path):
    num_polys = 0
    driver = ogr.GetDriverByName('ESRI Shapefile')
    hinputfile = driver.Open(input_path, 0)
    in_layer = hinputfile.GetLayer(0)
    for feature in in_layer:
        geom = feature.GetGeometryRef()
        if geom.GetGeometryType() == ogr.wkbPolygon:
            num_rings = geom.GetGeometryCount()
            num_points = 0
            num_polys += 1
            for i in range(num_rings):
                ring = geom.GetGeometryRef(i)
                num_points += ring.GetPointCount()
            print(num_polys, "\t", num_points)
        elif geom.GetGeometryType() == ogr.wkbMultiPolygon:
            for i in range(geom.GetGeometryCount()):
                g = geom.GetGeomteryRef(i)
                num_rings = g.GetGeometryCount()
                num_points = 0
                num_polys += 1
                for i in range(num_rings):
                    ring = g.GetGeometryRef(i)
                    num_points += ring.GetPointCount()
                print(num_polys, "\t", num_points)
        else:
            print("Unsupported Geometry Type!")
    hinputfile.Destroy()

if __name__ == "__main__":
    process_file(r'c:\GeoData\AZ.shp')

I only included the code for Polygon and Multi-Polygon types, but it is fairly easy to add the others.

Enjoy!