Python – Convert Polygon Feature Centroid to Points

Question

I would like to convert some polygon based shp files that have several polygon features into points for each feature that would essentially represent the centriod of each polygon feature. I know in ArcGIS world I could use Feature To Point tool but I would like to keep this in a script that could be run on PC's that don't have arcpy on them so I'm looking for an open source alternative to that. Is anyone aware of a library I could use for this along with some direction on how to leverage it accomplish this?

Answer 1

You can run an ogr2ogr command (e.g. from a OSGeo4w Shell). E.g. on a shapefile of countries:

cd path/to/shapefiles
ogr2ogr -sql "SELECT ST_Centroid(geometry), * FROM countries" -dialect sqlite countries_centroid.shp countries.shp

The new shapefile countries_centroid.shp should be similar to the input, but just contain one point per [Multi]Polygon.

@PEL also shows a good example with ST_PointOnSurface, which is simple to substitute in this command.

---

Something similar can be done in Python, if needed, but it may take a few lines of code more:

import os
from osgeo import ogr

ogr.UseExceptions()
os.chdir('path/to/shapefiles')

ds = ogr.Open('countries.shp')
ly = ds.ExecuteSQL('SELECT ST_Centroid(geometry), * FROM countries', dialect='sqlite')
drv = ogr.GetDriverByName('Esri shapefile')
ds2 = drv.CreateDataSource('countries_centroid.shp')
ds2.CopyLayer(ly, '')
ly = ds = ds2 = None  # save, close