import random as rd

stagelength = 4     # number of step between two verticille production
nbcycle = 7         # total number of verticille branches wanted 
minbranches = 2     # minimum number of branches when producing verticille
maxbranches = 5     # minimum number of branches when producing verticille
leafduration = 12   # life time of a leaf
leafmaturation = 6  # age at which a leaf is considered as old
angdev = 15         # authorized variation of angles to introduce some random
leafel = 60         # angle of leaf inclination from horizontal direction
radinc = 0.005      # increment of radius through time

def maxleafsize(s,maxs):      # maximum size of a leaf 
  return (s*0.5/float(maxs))+0.5

def leafsize(maxsize, t):     # size of a leaf according to its age
  if t < leafmaturation:
    return maxsize * ((t*0.5/float(leafmaturation))+0.5)
  else:
    return maxsize

def branch_angle(nc):         # branching angle according to position on trunk
    return 30+ 60*((nbcycle-nc)/float(nbcycle))

module A # represent trunk apical meristem
module B # represent apical meristem of lateral branches
module L # whorl of leaf on trunk
module K # whorl of leaf on lateral branches
module I # Internode

Axiom: SetWidth(0.1)I(0.5,0.1)A(-stagelength+1)

derivation length: nbcycle*stagelength
production:

A(t) :
  nproduce I(1,0.1)/(90)L(maxleafsize(t%stagelength,stagelength)*1.5,0) # produce a metamer
  if (t % stagelength) == 0:
      cyclenb = t // stagelength
      # produce a verticille of branches.
      nbaxe = rd.randint(3,7)
      for i in range(nbaxe):
        nproduce [/(360*i/nbaxe)&(branch_angle(cyclenb))B]
  produce A(t+1)
  # produce apex with older age.

B --> I(0.5,0.1)K(0)B

K(t) :
  # ages the branches leaves. If too old, removes
  if t <= leafduration :
    produce K(t+1)
  else:
    produce

L(maxsize,t) :
  # ages the trunk leaves. If too old, removes
  if t <= leafduration :
    produce L(maxsize,t+1)
  else:
    produce *  

I(s,r) --> I(s,r+radinc)


homomorphism:

I(a,r) --> SetWidth(r)F(a, r-radinc)
L(maxsize,t) --> [&(120);(2)~l(leafsize(maxsize,t))][^(120);(2)~l(leafsize(maxsize,t))]

K(t) --> [+(rd.uniform(90-angdev, 90+angdev))&(rd.uniform(leafel-angdev, leafel+angdev));(2)~l(leafsize(1,t))][-(rd.uniform(90-angdev, 90+angdev))&(rd.uniform(leafel-angdev, leafel+angdev));(2)~l(leafsize(1,t))]

endlsystem

###### INITIALISATION ######

def __initialiseContext__(context):
	context.animation_timestep = 0.153
	context.options.setSelection("Module declaration",1)
__description__ = """Massart architectural model."""
__authors__ = """E. COSTES, F. Boudon"""
