# How to create a dense pointcloud of trees in ArcGIS?

I'm having a large dataset of forest polygons and would like to scatter individual tree locations over these areas. Here is how I have done this so far: In "Create Random Points" number of points has been set to 999999 and the minimum allowed distance (linear) to 5 meters.

Why am I running the "create random points" 4 times in that model? The problem is that when running it once in ArcGIS, the points created are way too sparse. Running it four times and then merging all positions that are closer together than 1.4 meters (which essentially comes out to ~5 meters as tests have shown) gives me a much denser point cloud.

Any ideas how to change this approach into something faster with improved performance?

• Is this purely for visualisation purposes? Or do you need auto-generated point features to run some kind of hypothetical analysis? – RobinHorner Aug 29 '18 at 5:02
• It's for visualization (placement of 3D tree models) – Brief Aug 29 '18 at 5:38

You dont mention how fast your model executes. I tried code below for three polygons, generating about 60000 random trees in one minute. Make sure not to input impossible combination of tree density and minimum tree distance or it will run forever.

``````import arcpy
from random import randint
arcpy.env.overwriteOutput = True

#Change these five lines
arcpy.env.workspace = r'C:\Folder\Default.gdb'
stands = r'Stands'
output_points = r'Randomtrees' #Will be created in script
trees_per_areaunit = 0.1 #Desired tree density
min_treedistance = 2 #Min distance between trees

arcpy.CreateFeatureclass_management(out_path=arcpy.env.workspace, out_name=output_points, geometry_type='POINT', spatial_reference=stands)

def giverandomcoord(xmin,xmax,ymin,ymax):
return (randint(xmin,xmax),randint(ymin,ymax))

def givenogoarea(x, y, sigma):
#https://stackoverflow.com/questions/39862709/generate-coordinates-in-grid-that-lie-within-a-circle
neighborhood = []
X = int(sigma)
for i in range(-X, X + 1):
Y = int(pow(sigma * sigma - i * i, 1/2))
for j in range(-Y, Y + 1):
neighborhood.append((x + i, y + j))
return neighborhood

with arcpy.da.SearchCursor(stands,['SHAPE@','OID@']) as cursor:
for row in cursor:
coords = []
nogo = set()
ext = row.extent
treegoal = int(ext.polygon.area*trees_per_areaunit)
xmin, xmax, ymin, ymax = int(ext.lowerLeft.X), int(ext.upperRight.X), int(ext.lowerLeft.Y), int(ext.upperRight.Y)
while len(coords)<treegoal:
point = giverandomcoord(xmin, xmax, ymin, ymax)
if point not in nogo:
coords.append(point)
nogo.update(givenogoarea(point, point, min_treedistance))
pointlist = [arcpy.PointGeometry(arcpy.Point(*c)) for c in coords]
arcpy.CopyFeatures_management(in_features=pointlist, out_feature_class=r'in_memory\temppoints')
arcpy.MakeFeatureLayer_management(in_features=r'in_memory\temppoints', out_layer='templyrpoints')
arcpy.MakeFeatureLayer_management(in_features=stands, out_layer=r'templyrstand', 