Generating a Grid of Points on Selected Polygons

[MatthewPurdy]

Just wondering if anyone out there has any ideas on how to generate a systematic grid of points that are a specified distance apart via python. I'm creating plots for field measurements and for each polygon that represents a forested stand, I'm creating a regularly spaced grid for sampling purposes. I'm aware of third party tools (i.e.Hawths) but I'd like to do it in python. I've looked at the 'create random points' and 'fishnet' tools but can't really think of how to impliment them in this manner. I'm using v10. Let me know if you need some clarifications.

thanks.

[DarrenWiens2]

1.) Polygon to Raster tool, with cell size for your grid.
2.) Raster to Point tool, to get your point grid.

[ChrisSnyder]

Question: Do you want the points like this

#aligned like a square (not equadistant)
+   +   +   +

+   +   +   +

+   +   +   +

+   +   +   +

or like this:

#aligned like a tringle/diamond (equadistant)
+   +   +   +
  
   +   +   +   +

+   +   +   +

   +   +   +   +

If it's the former, just use the fishnet tool (and convert the polygon centroids to points using the FeatureToPoint tool), or as Darren suggests, use a raster/vector hybrid approach  If it's the latter (which is what we use BTW to lay out our inventory plots), I think I have some code that I can borrow from someone. Let me know, and I'll try to post it (if they let me).

[JosefRundström1]

Oh! I would loooove to see that small piece of code!

[SHoldini]

Question: Do you want the points like this

#aligned like a square (not equadistant)
+   +   +   +

+   +   +   +

+   +   +   +

+   +   +   +

or like this:

#aligned like a tringle/diamond (equadistant)
+   +   +   +
  
   +   +   +   +

+   +   +   +

   +   +   +   +

If it's the former, just use the fishnet tool (and convert the polygon centroids to points using the FeatureToPoint tool), or as Darren suggests, use a raster/vector hybrid approach  If it's the latter (which is what we use BTW to lay out our inventory plots), I think I have some code that I can borrow from someone. Let me know, and I'll try to post it (if they let me).

Hi, I would also be interested in seeing the code for the approach to this, thanks:)

[ChrisSnyder]

Late is better than never....

This code generates equadistant points that more than covers the extent out the 'containerFC'. Further code could be written to limit the points to those that are contained in the containerFC... either in the insert cursor using the .contains() geometry method, or at the very end, a SelectByLocation and then export to a new FC.

I put in a 'useRandomShift' flag so that the points don't necessarily have to be tied to the containerFC's extent.

#Chris Snyder, WADNR, 20120126
import arcpy, math, random
containerFC = r"C:\csny490\test.shp"
outputFC = r"C:\csny490\out.shp"
spacing = 100
useRandomShift = True
arcpy.env.overwriteOutput = True
dsc = arcpy.Describe(containerFC)
dscExtent = dsc.extent
xMin = dscExtent.xMin
yMin = dscExtent.yMin
xMax = dscExtent.xMax
yMax = dscExtent.yMax
randomXshift = random.uniform(0, spacing)
randomYshift = random.uniform(0, spacing)
yOffset = math.sqrt((spacing **2) - ((spacing / 2) **2))
pntDict = {}
i = 0
rowNumber = 1
xCoord = xMin - spacing
yCoord = yMin - spacing
while xCoord < xMax + spacing or yCoord < yMax + spacing:
    i = i + 1
    if xCoord > xMax + spacing:
        rowNumber = rowNumber + 1
        xCoord = xMin - spacing
        if bool(rowNumber & 1) == False:
            xCoord = xCoord + spacing / 2
        yCoord = yCoord + yOffset
    elif i == 1:
        xCoord = xCoord
    else:
        xCoord = xCoord + spacing
    if useRandomShift == True:
        pntDict = (xCoord + randomXshift, yCoord + randomYshift)
    else:
        pntDict = (xCoord, yCoord)
arcpy.CreateFeatureclass_management(outputFC[:-len(outputFC.split("\\")[-1]) - 1], outputFC.split("\\")[-1], "POINT", "", "", "", dsc.spatialreference, "", "", "", "")
shapeFieldName = arcpy.Describe(outputFC).shapeFieldName
insertRows = arcpy.InsertCursor(outputFC)
for pnt in pntDict:
    insertRow = insertRows.newRow()
    pointObj = arcpy.CreateObject("Point")
    pointObj.X = pntDict[pnt][0]
    pointObj.Y = pntDict[pnt][1]
    insertRow.setValue(shapeFieldName, pointObj)
    insertRows.insertRow(insertRow)
del pointObj, insertRow, insertRows