8

I have a set of polygons I have calculated the centroids of, but am having trouble finding a way to output these values to a point shapefile. I dabbled with arcpy.MakeXYEvents_management, but it only creates a temporary layer file.

4
  • 7
    you can export that XY event table to a new feature class.
    – crmackey
    Dec 16, 2015 at 20:50
  • Wow! A lot of solutions. I am in. This gis.stackexchange.com/questions/147790/… will create the centre of polygon's largest inscribed circle and transfer all visible fields into output point feature class. A bit slow though...
    – FelixIP
    Dec 16, 2015 at 23:32
  • I'll have to remember that @FelixIP, that code is quite nice and I can see situations where it would come in handy. Dec 17, 2015 at 1:07
  • @MichaelMiles-Stimson thanks. Definitely helped here gis.stackexchange.com/questions/166230/… although it won't work on donuts and multipart polygons
    – FelixIP
    Dec 17, 2015 at 1:19

4 Answers 4

11

Here is an arcpy soloution:

import os, sys, arcpy

InFC = sys.argv[1]
OutFC = sys.argv[2]

#split the output into directory and name
Folder = os.path.dirname(OutFC)
Name   = os.path.basename(OutFC)

# Get the existins spatial reference as it's going to match
desc = arcpy.Describe(InFC)
SR   = desc.spatialReference

# create or append
if not arcpy.Exists(OutFC):
    arcpy.CreateFeatureclass_management(Folder,Name,"POINT",spatial_reference = SR)

# open up the 'writer'
with arcpy.da.InsertCursor(OutFC,"SHAPE@XY") as iCur:
    # open up the 'reader'
    with arcpy.da.SearchCursor(InFC,"SHAPE@") as sCur:
        # loop through each polygon in the InFC
        for sRow in sCur:
            cent = sRow[0].centroid          # get the centroid
            iCur.insertRow([(cent.X,cent.Y)])# write it to the new feature class

This creates a new feature class, opens an insert cursor, loops through each geometry in the InFC and writes the centroid property to the output.. note that multi-part geometries will not have a centroid per part but rather a single centroid for the whole geometry.. it's not much more difficult to loop through the parts - let's keep this example fairly simple though to show the basics.

1
  • +1 for beating me to the punch with a nice concise example!
    – crmackey
    Dec 16, 2015 at 21:41
8

This should do the trick:

import arcpy
import os

def polysToPoints(in_polys, out_points):
    """converts polygons to centroids

    in_polys -- input polygons
    out_points -- output points
    """
    ws, name = os.path.split(out_points)
    sr = arcpy.Describe(in_polys).spatialReference
    arcpy.management.CreateFeatureclass(ws, name, 'POINT', template=in_polys, spatial_reference=sr)

    # populate records
    fields = ['SHAPE@'] + [f.name for f in arcpy.ListFields(out_points) if f.type not in ('OID', 'Geometry')]
    with arcpy.da.InsertCursor(out_points, fields) as irows:
        with arcpy.da.SearchCursor(in_polys, fields) as rows:
            for r in rows:
                irows.insertRow((arcpy.PointGeometry(r[0].centroid),) + r[1:])
    return out_points

if __name__ == '__main__':

    polys = r'C:\path_to_your\polys.shp'
    points = r'C:\path_to_your\points.shp'
    polysToPoints(polys, points)
3
  • Haha wow, @Michael Miles-Stimson posted something very similar right before me...The only thing mine does differently is adds all the field values with each record.
    – crmackey
    Dec 16, 2015 at 21:40
  • 1
    I like that, it's good to have all the fields copied... +1 from me. How about multi-part polygons/polylines? They aren't common, except for holes, but appear often enough to cause problems. Dec 17, 2015 at 1:06
  • Yeah, my example does not account for that. It would be easy enough though to iterate through each part (which returns an array) build a polygon object from that array and generate a centroid from that and append to a geometry for a multipoint feature class.
    – crmackey
    Dec 17, 2015 at 20:52
4

shapelib library is ideal for this task.
No depedecies, just a simple module to create shapefiles:

import shapefile

coord_list = ((20,40), (0,0), )

w = shapefile.Writer(shapefile.POINT)
w.field('ID_FIELD','C','40')
for id,(x,y) in enumerate(coord_list):
    w.record(ID_FIELD=id)
    w.point(x,y)
w.save("file\to\save")

it's very simple you might want to add the projection file afterwards (manually).

Check the manual on additional info

1
  • Wow lot's of responses, I'll give all these great suggestions a shot tomorrow and update as to which ran optimally for this project. Thanks all!
    – geodranic
    Dec 17, 2015 at 1:29
1

I'm still somewhat of a an arcpy beginner so perhaps this is not as sophisticated, but I created a script tool for this and it works for me. In your toolbox, add a new script tool. Make the python script below the source. Set the tools 1st parameter as feature layer input and the 2nd parameter as feature layer output (make sure you set it to an output parameter). Note, this tool will return points inside the input layer rather than the true centroid.

import arcpy
import os

inputFC = arcpy.GetParameterAsText(0)
outputPath = arcpy.GetParameterAsText(1)
scratchPath = arcpy.env.scratchGDB
inputFCSpatialReference = arcpy.Describe(inputFC).spatialReference

##copy to new fc in scratchGBD and make new layer
FCCopy = arcpy.CopyFeatures_management(inputFC,os.path.join(scratchPath,"FCCopy"))
FCLyr = arcpy.MakeFeatureLayer_management(FCCopy,"FCLyr")
arcpy.AddMessage("FC LAYER CREATED")

##add and calculate X and Y fields
arcpy.AddField_management(FCLyr,"Xcoor","DOUBLE")
arcpy.CalculateField_management(FCLyr,"Xcoor","!SHAPE!.centroid.X","PYTHON")
arcpy.AddField_management(FCLyr,"Ycoor","DOUBLE")
arcpy.CalculateField_management(FCLyr,"Ycoor","!SHAPE!.centroid.Y","PYTHON")
arcpy.AddMessage("X,Y FIELDS CREATED")

##create point layer from X Y fields
FCTable = arcpy.MakeTableView_management(FCLyr,"FCTable")
pointsLyr = arcpy.MakeXYEventLayer_management(FCTable,"Xcoor","Ycoor","pointsLyr",inputFCSpatialReference)
arcpy.AddMessage("TABLE LAYER CREATED")

##export points to feature class
arcpy.CopyFeatures_management(pointsLyr,outputPath)
arcpy.AddMessage("POINT LAYER CREATED")

##delete data from scratchGBD
arcpy.Delete_management(FCCopy)

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.