# Create a 3D tree based on 3D polygons with Arcpy

Blog Post created by xander_bakker on Jul 10, 2016

At this point I suppose you already read my posts on Working with 3D and M-aware geometries in Arcpy and Creating a 3D tree with Arcpy for 3D Analysis. The idea is to create a 3D tree using polygons. Something like the one below (which is only partially shown):

My first instinct was to take the code used to create the multipart polylines and divide these in parts and connect the parts from two sequential rings into a polygon. The first method I tried was to take the buffer outline and use the Polyline::positionAlongLine (value, {use_percentage}) method to pick the points on the lines. When using the same value these point should align on each ring. At least, that was what I thought. However, that was no success, leading to something very artistic, but little useful.

A very easy way of doing something similar is using the PointGeometry::pointFromAngleAndDistance (angle, distance, {method}). This yields a point based on an input point, an angle and a distance. Sweet!

The second challenge was changing the initial though of creating a multipart geometry, one for the trunk and one for the crown. This is no problem when you create a multipart polyline, but when you create a multipart polygon, it will apply some 2D topology rules to create the multipart polygon. This results in the loss of many polygon parts:

The third aspect you should bear in mind is that vertical polygons will not work in Arcpy. Hence the idea: (Better) support for Multipatches in Arcpy . This can be easily avoided by letting the trunk decrease in diameter as it extends upwards and by avoiding an even number of crown_circles.

When you use an even number, you get this:

... a tree with a some kind of advanced ventilation implemented.

So let's see what we can do to create the tree using polygons (not multipart) and throw in some python code:

```import arcpy
import math

def main():
import os

# output
fc_out = r'C:\GeoNet\_WorkingWith3D\data.gdb\TreePolygon3D_v01'
arcpy.env.overwriteOutput = True

# the tree point
sr = arcpy.SpatialReference(3857)
x = -13046027.327
y = 4036381.962
elevation = 365

# dimensions of tree
tree_height = 30
crown_diameter = 16
trunk_diameter = 0.7
trunk_height = 15
crown_height = tree_height - trunk_height

# configuration of number of points and lines (detail of tree)
crown_lines = 15  # number of lines for crown (avoid even numbers!)
trunk_lines = 25 # number of lines for trunk
pnts_on_circle = 18  # a point every 20 degrees on circle

# create output featureclass
ws, fc_name = os.path.split(fc_out)
arcpy.CreateFeatureclass_management(ws, fc_name, "POLYGON", None, None, "ENABLED", sr)

# Addional field for tree part to store Trunk or Crown
fld_tree_part = "TreePart"
arcpy.AddField_management(fc_out, fld_tree_part, "TEXT", None, None, 10)

# insert cursor to write feature
flds_out = ('SHAPE@', fld_tree_part)
with arcpy.da.InsertCursor(fc_out, flds_out) as curs_out:

# list to store multipart feature
lst_feats = []
pnt = arcpy.Point(x, y, elevation)
pntg = arcpy.PointGeometry(pnt, sr, True)

# define z values
z0 = elevation
z1 = z0 + trunk_height
z2 = z1 + crown_height

# create list of lists with points to create trunk polygons
lst_polygons = createTrunk3Dv02(pntg, trunk_lines, pnts_on_circle, z0, z1 +1, trunk_diameter, sr)
for lst_pol in lst_polygons:
pol_trunk = arcpy.Polygon(arcpy.Array(lst_pol), sr, True, False)
row_out = (pol_trunk, "Trunk", )
curs_out.insertRow(row_out)

# create list of lists with points to create crown polygons
lst_polygons = createCrown3Dv02(pntg, crown_lines, pnts_on_circle, z1, z2, crown_diameter, sr)
for lst_pol in lst_polygons:
pol_trunk = arcpy.Polygon(arcpy.Array(lst_pol), sr, True, False)
row_out = (pol_trunk, "Crown", )
curs_out.insertRow(row_out)

def createTrunk3Dv02(pntg, trunk_lines, pnts_on_circle, z0, z1, trunk_diameter, sr):
'''Create a list of lists where each inner list contains the 5 points
to create a 3D polygon as part of the trunk'''

lst_polygons = []
for a in range(pnts_on_circle):
# convert point on circle to angle 0 - 360 degrees
angle1 = float(a) / float(pnts_on_circle) * 360
angle2 = float(a + 1) / float(pnts_on_circle) * 360

# create 2 points using first and second angle and trunk diameter
pntga1d1 = pntg.pointFromAngleAndDistance(angle1, trunk_diameter, 'PLANAR')
pntga2d1 = pntg.pointFromAngleAndDistance(angle2, trunk_diameter, 'PLANAR')

# create 2 points using first and second angle and 90% of trunk diameter
pntga1d2 = pntg.pointFromAngleAndDistance(angle1, trunk_diameter * 0.9, 'PLANAR')
pntga2d2 = pntg.pointFromAngleAndDistance(angle2, trunk_diameter * 0.9, 'PLANAR')

# make points z-aware
pnta1d1 = arcpy.Point(pntga1d1.firstPoint.X, pntga1d1.firstPoint.Y, z0)
pnta2d1 = arcpy.Point(pntga2d1.firstPoint.X, pntga2d1.firstPoint.Y, z0)
pnta1d2 = arcpy.Point(pntga1d2.firstPoint.X, pntga1d2.firstPoint.Y, z1)
pnta2d2 = arcpy.Point(pntga2d2.firstPoint.X, pntga2d2.firstPoint.Y, z1)

# add 4 points + 1st point as list to polygon list
lst_polygons.append([pnta1d1, pnta2d1, pnta2d2, pnta1d2, pnta1d1])

# return list with lists of points to create polygons
return lst_polygons

def createCrown3Dv02(pntg, crown_lines, pnts_on_circle, z1, z2, crown_diameter, sr):
'''Create a list of lists where each inner list contains the 5 points
to create a 3D polygon as part of the crown'''

lst_polygons = []
for n in range(crown_lines):
# take pairs of rings (n and n+1)
i1 = n
i2 = n + 1

# determine the elevation for each of the two rings
zi1 = float(i1) / float(crown_lines) * (z2  - z1) + z1
zi2 = float(i2) / float(crown_lines) * (z2  - z1) + z1

# convert i1 and i2 to a value expressed in pi
p1 = float(i1) / (crown_lines + 1) * math.pi
p2 = float(i2) / (crown_lines + 1) * math.pi
# calculate the diameter based on this value expressed in pi
diam1 = math.sin(p1) * float(crown_diameter)
diam2 = math.sin(p2) * float(crown_diameter)

for a in range(pnts_on_circle):
# determine the start and end angle for each pair of points on circle
angle1 = float(a) / float(pnts_on_circle) * 360
angle2 = float(a + 1) / float(pnts_on_circle) * 360

# create 4 points, two on each ring
pntga1d1 = pntg.pointFromAngleAndDistance(angle1, diam1 / 2.0, 'PLANAR')
pntga2d1 = pntg.pointFromAngleAndDistance(angle2, diam1 / 2.0, 'PLANAR')
pntga1d2 = pntg.pointFromAngleAndDistance(angle1, diam2 / 2.0, 'PLANAR')
pntga2d2 = pntg.pointFromAngleAndDistance(angle2, diam2 / 2.0, 'PLANAR')

# make points z-aware
pnta1d1 = arcpy.Point(pntga1d1.firstPoint.X, pntga1d1.firstPoint.Y, zi1)
pnta2d1 = arcpy.Point(pntga2d1.firstPoint.X, pntga2d1.firstPoint.Y, zi1)
pnta1d2 = arcpy.Point(pntga1d2.firstPoint.X, pntga1d2.firstPoint.Y, zi2)
pnta2d2 = arcpy.Point(pntga2d2.firstPoint.X, pntga2d2.firstPoint.Y, zi2)

# add 4 points + 1st point as list to polygon list
lst_polygons.append([pnta1d1, pnta2d1, pnta2d2, pnta1d2, pnta1d1])

# return list with lists of points to create polygons
return lst_polygons

if __name__ == '__main__':
main()
```

If you have any questions about the code, just add a comment and I will get back to you.

Also have a look at the 3D City Base Layers - Add Trees topic of the 3D Cities section, there is a tool to create trees: