1- # Forest
21
3- Alejandro Morales and Ana Ernst
42
5- Centre for Crop Systems Analysis - Wageningen University
3+ ```` @example Forest
4+ #= Forest
5+
6+ Alejandro Morales
67
7- > ## TL;DR
8- > Similar in functionality to [ Tree] ( https://virtualplantlab.com/dev/tutorials/from_tree_forest/tree/ ) tutorial with separate graphs for each tree
9- > - Modify tree parameters for each tree
10- > - Multithreaded simulation (grow trees in parallel)
11- > - Scene customization (e.g., add soil)
12- > - Export Scenes
8+ Centre for Crop Systems Analysis - Wageningen University
139
1410In this example we extend the tree example into a forest, where
1511each tree is described by a separate graph object and parameters driving the
1612growth of these trees vary across individuals following a predefined distribution.
13+ This tutorial requires using the Distributions.jl package:
14+
1715The data types, rendering methods and growth rules are the same as in the binary
1816tree example:
1917
20- ```` julia
18+ =#
2119using VirtualPlantLab
2220using Distributions, Plots, ColorTypes
2321import GLMakie
@@ -28,8 +26,8 @@ module TreeTypes
2826 struct Meristem <: VirtualPlantLab.Node end
2927 # Bud
3028 struct Bud <: VirtualPlantLab.Node end
31- # TreeNode
32- struct TreeNode <: VirtualPlantLab.Node end
29+ # Node
30+ struct Node <: VirtualPlantLab.Node end
3331 # BudNode
3432 struct BudNode <: VirtualPlantLab.Node end
3533 # Internode (needs to be mutable to allow for changes over time)
@@ -46,29 +44,23 @@ module TreeTypes
4644 growth::Float64 = 0.1
4745 budbreak::Float64 = 0.25
4846 phyllotaxis::Float64 = 140.0
49- leaf_angle:: Float64 = 45 .0
50- branch_angle:: Float64 = 30 .0
47+ leaf_angle::Float64 = 30 .0
48+ branch_angle::Float64 = 45 .0
5149 end
5250end
5351
5452import .TreeTypes
55- ````
5653
57- Create geometry and color for the * internodes* :
58-
59- ```` julia
54+ # Create geometry + color for the internodes
6055function VirtualPlantLab.feed!(turtle::Turtle, i::TreeTypes.Internode, data)
6156 # Rotate turtle around the head to implement elliptical phyllotaxis
6257 rh!(turtle, data.phyllotaxis)
6358 HollowCylinder!(turtle, length = i.length, height = i.length/15, width = i.length/15,
6459 move = true, colors = RGB(0.5,0.4,0.0))
6560 return nothing
6661end
67- ````
68-
69- Create geometry and color for the * leaves* :
7062
71- ```` julia
63+ # Create geometry + color for the leaves
7264function VirtualPlantLab.feed!(turtle::Turtle, l::TreeTypes.Leaf, data)
7365 # Rotate turtle around the arm for insertion angle
7466 ra!(turtle, -data.leaf_angle)
@@ -79,21 +71,16 @@ function VirtualPlantLab.feed!(turtle::Turtle, l::TreeTypes.Leaf, data)
7971 ra!(turtle, data.leaf_angle)
8072 return nothing
8173end
82- ````
83-
84- Insertion angle for the bud nodes:
8574
86- ```` julia
75+ # Insertion angle for the bud nodes
8776function VirtualPlantLab.feed!(turtle::Turtle, b::TreeTypes.BudNode, data)
8877 # Rotate turtle around the arm for insertion angle
8978 ra!(turtle, -data.branch_angle)
9079end
91- ````
9280
93- Rules for meristem and branches:
9481
95- ```` julia
96- meristem_rule = Rule (TreeTypes. Meristem, rhs = mer -> TreeTypes. TreeNode () +
82+ # Rules
83+ meristem_rule = Rule(TreeTypes.Meristem, rhs = mer -> TreeTypes.Node () +
9784 (TreeTypes.Bud(), TreeTypes.Leaf()) +
9885 TreeTypes.Internode() + TreeTypes.Meristem())
9986
@@ -129,7 +116,7 @@ of each tree and rotates it.
129116(ii) Wrap the axiom, rules and the creation of the graph into a function that
130117takes the required parameters as inputs.
131118
132- ```` julia
119+ ```` @example Forest
133120function create_tree(origin, growth, budbreak, orientation)
134121 axiom = T(origin) + RH(orientation) + TreeTypes.Internode() + TreeTypes.Meristem()
135122 tree = Graph(axiom = axiom, rules = (meristem_rule, branch_rule),
141128The code for elongating the internodes to simulate growth remains the same as for
142129the binary tree example
143130
144- ```` julia
131+ ```` @example Forest
145132getInternode = Query(TreeTypes.Internode)
146133
147134function elongate!(tree, query)
@@ -169,43 +156,50 @@ of 2 meters. First we generate the original positions of the trees. For the
169156position we just need to pass a ` Vec ` object with the x, y, and z coordinates of
170157the location of each TreeTypes. The code below will generate a matrix with the coordinates:
171158
172- ```` julia
159+ ```` @example Forest
173160origins = [Vec(i,j,0) for i = 1:2.0:20.0, j = 1:2.0:20.0]
174161````
175162
176163We may assume that the initial orientation is uniformly distributed between 0 and 360 degrees:
177164
178- ```` julia
165+ ```` @example Forest
179166orientations = [rand()*360.0 for i = 1:2.0:20.0, j = 1:2.0:20.0]
180167````
181168
182- For the ` growth ` and ` budbreak ` parameters we will be assumed to follow a
169+ For the ` growth ` and ` budbreak ` parameters we will assumed that they follow a
183170LogNormal and Beta distribution, respectively. We can generate random
184171values from these distributions using the ` Distributions ` package. For the
185172relative growth rate:
186173
187- ```` julia
188- growths = rand (LogNormal (- 2 , 0.2 ), 10 , 10 )
174+ ```` @example Forest
175+ growths = rand(LogNormal(-2, 0.3 ), 10, 10)
189176histogram(vec(growths))
177+ savefig("growths.png") ## hide
190178````
191179
180+ ![ ] ( growths.png )
181+
192182And for the budbreak parameter:
193183
194- ```` julia
184+ ```` @example Forest
195185budbreaks = rand(Beta(2.0, 10), 10, 10)
196186histogram(vec(budbreaks))
187+ savefig("budbreaks.png") ## hide
197188````
198189
190+ ![ ] ( budbreaks.png )
191+
199192Now we can create our forest by calling the ` create_tree ` function we defined earlier
200193with the correct inputs per tree:
201194
202- ```` julia
203- forest = vec (create_tree .(origins, growths, budbreaks, orientations))
195+ ```` @example Forest
196+ forest = vec(create_tree.(origins, growths, budbreaks, orientations));
197+ nothing #hide
204198````
205199
206200By vectorizing ` create_tree() ` over the different arrays, we end up with an array
207201of trees. Each tree is a different Graph, with its own nodes, rewriting rules
208- and variables. This avoids having to create large graphs to include all the
202+ and variables. This avoids having to create a large graphs to include all the
209203plants in a simulation. Below we will run a simulation, first using a sequential
210204approach (i.e. using one core) and then using multiple cores in our computers (please check
211205https://docs.julialang.org/en/v1/manual/multi-threading/ if the different cores are not being used
@@ -216,108 +210,126 @@ as you may need to change some settings in your computer).
216210We can simulate the growth of each tree by applying the method ` simulate ` to each
217211tree, creating a new version of the forest (the code below is an array comprehension)
218212
219- ```` julia
220- newforest = [simulate (tree, getInternode, 2 ) for tree in forest]
213+ ```` @example Forest
214+ newforest = [simulate(tree, getInternode, 2) for tree in forest];
215+ nothing #hide
221216````
222217
223218And we can render the forest with the function ` render ` as in the binary tree
224219example but passing the whole forest at once
225220
226- ```` julia
227- render (Scene (newforest))
221+ ```` @example Forest
222+ pl = render(Scene(newforest))
223+ GLMakie.save("newforest1.png", pl) ## hide
228224````
229225
226+ ![ ] ( newforest1.png )
227+
230228If we iterate 4 more iterations we will start seeing the different individuals
231229diverging in size due to the differences in growth rates
232230
233- ```` julia
231+ ```` @example Forest
234232newforest = [simulate(tree, getInternode, 15) for tree in newforest];
235- render (Scene (newforest))
233+ pl = render(Scene(newforest))
234+ GLMakie.save("newforest2.png", pl) ## hide
236235````
237236
237+ ![ ] ( newforest2.png )
238+
238239## Multithreaded simulation
239240
240241In the previous section, the simulation of growth was done sequentially, one tree
241242after another (since the growth of a tree only depends on its own parameters). However,
242243this can also be executed in multiple threads. In this case we use an explicit loop
243244and execute the iterations of the loop in multiple threads using the macro ` @threads ` .
244- Note that the rendering function can also be run in parallel (i.e. the geometry will be
245+ Note that the rendering function can also be ran in parallel (i.e. the geometry will be
245246generated separately for each plant and the merge together):
246247
247- ```` julia
248+ ```` @example Forest
248249using Base.Threads
249250newforest = deepcopy(forest)
250- @threads for i in eachindex (forest)
251+ @threads for i in 1:length (forest)
251252 newforest[i] = simulate(forest[i], getInternode, 6)
252253end
253- render (Scene (newforest), parallel = true )
254+ pl = render(Scene(newforest, parallel = true))
255+ GLMakie.save("newforest3.png", pl) ## hide
254256````
255257
258+ ![ ] ( newforest3.png )
259+
256260An alternative way to perform the simulation is to have an outer loop for each timestep and an internal loop over the different trees. Although this approach is not required for this simple model, most FSP models will probably need such a scheme as growth of each individual plant will depend on competition for resources with neighbouring plants. In this case, this approach would look as follows:
257261
258- ```` julia
262+ ```` @example Forest
259263newforest = deepcopy(forest)
260264for step in 1:15
261- @threads for i in eachindex (newforest)
265+ @threads for i in 1:length (newforest)
262266 newforest[i] = simulate(newforest[i], getInternode, 1)
263267 end
264268end
265- render (Scene (newforest), parallel = true )
269+ pl = render(Scene(newforest, parallel = true))
270+ GLMakie.save("newforest4.png", pl) ## hide
266271````
267272
273+ ![ ] ( newforest4.png )
274+
268275# Customizing the scene
269276
270- Here we are going to customize the scene of our simulation by adding a horizontal tile representing soil and
277+ Here we are going to customize the scene of our simulation by adding a horizontal tile represting soil and
271278tweaking the 3D representation. When we want to combine plants generated from graphs with any other
272279geometric element it is best to combine all these geometries in a ` GLScene ` object. We can start the scene
273280with the ` newforest ` generated in the above:
274281
275- ```` julia
276- scene = Scene (newforest)
282+ ```` @example Forest
283+ scene = Scene(newforest);
284+ nothing #hide
277285````
278286
279287We can create the soil tile directly, without having to create a graph. The simplest approach is two use
280288a special constructor ` Rectangle ` where one species a corner of the rectangle and two vectors defining the
281289two sides of the vectors. Both the sides and the corner need to be specified with ` Vec ` just like in the
282290above when we determined the origin of each plant. VPL offers some shortcuts: ` O() ` returns the origin
283- (` Vec(0.0, 0.0, 0.0) ` ), whereas ` X ` , ` Y ` and ` Z ` returns the corresponding axes, and you can scale them by
291+ (` Vec(0.0, 0.0, 0.0) ` ), whereas ` X ` , ` Y ` and ` Z ` returns the corresponding axes and you can scale them by
284292passing the desired length as input. Below, a rectangle is created on the XY plane with the origin as a
285293corner and each side being 11 units long:
286294
287- ```` julia
295+ ```` @example Forest
288296soil = Rectangle(length = 21.0, width = 21.0)
289297rotatey!(soil, pi/2)
290298VirtualPlantLab.translate!(soil, Vec(0.0, 10.5, 0.0))
291299````
292300
293301We can now add the ` soil ` to the ` scene ` object with the ` add! ` function.
294302
295- ```` julia
303+ ```` @example Forest
296304VirtualPlantLab.add!(scene, mesh = soil, colors = RGB(1,1,0))
297305````
298306
299307We can now render the scene that combines the random forest of binary trees and a yellow soil. Notice that
300308in all previous figures, a coordinate system with grids was being depicted. This is helpful for debugging
301- your code but also to help set up the scene (e.g. if you are not sure how big the soil tile should be).
302- However , it may be distracting for the visualization. It turns out that we can turn that off with
303- ` show_axes = false` :
309+ your code but also to help setup the scene (e.g. if you are not sure how big the soil tile should be).
310+ Howver , it may be distracting for the visualization. It turns out that we can turn that off with
311+ ` axes = false` :
304312
305- ```` julia
306- render (scene, axes = false )
313+ ```` @example Forest
314+ pl = render(scene, axes = false)
315+ GLMakie.save("newforest5.png", pl) ## hide
307316````
308317
318+ ![ ] ( newforest5.png )
319+
309320We may also want to save a screenshot of the scene. For this, we need to store the output of the ` render ` function.
310321We can then resize the window rendering the scene, move around, zoom, etc. When we have a perspective that we like,
311322we can run the ` save_scene ` function on the object returned from ` render ` . The argument ` resolution ` can be adjusted in both
312323` render ` to increase the number of pixels in the final image. A helper function ` calculate_resolution ` is provided to
313324compute the resolution from a physical width and height in cm and a dpi (e.g., useful for publications and posters):
314325
315- ```` julia
326+ ```` @example Forest
316327res = calculate_resolution(width = 16.0, height = 16.0, dpi = 1_000)
317- output = render (scene, axes = false , resolution = res)
328+ output = render(scene, axes = false, size = res)
318329export_scene(scene = output, filename = "nice_trees.png")
319330````
320331
321332---
322333
323334* This page was generated using [ Literate.jl] ( https://github.com/fredrikekre/Literate.jl ) .*
335+
0 commit comments