Sunday, January 15, 2012

Of Trees and Turtles

An L-system is a fascinating crossbreed between formal language theory, botany and computer graphics. By studying how it works, it's really not hard to imagine that nature, however it does it, must be doing something at least remotely similar. It's really easy to create an L-system rewrite engine in Python:

rules = {'A': '[>FA][^FA]'}
axiom = 'FA'
n_iters = 2

seq = [''.join([rules.get(c,c) 
                for c in (locals()['_[1]'][-1] if locals()['_[1]'] else axiom)]) for i in range(n_iters)][-1]

print seq # F[>F[>FA][^FA]][^F[>FA][^FA]]

A simple one-liner, really? Wow, isn't Python nice and easy? But wait.. it's wrong! I was kidding, because this is an absolutely ridiculous way of doing it, using an obscure Python feature allowing a list to refer to itself.. really not recommended! This is a less spectacular, but much saner way of doing it:

seq = axiom
for i in range(n_iters):
    rewrite = ''
    for c in seq:
        rewrite += rules.get(c, c)
    seq = rewrite
    
print seq # F[>F[>FA][^FA]][^F[>FA][^FA]]

We've "grown" a string, but now how about doing something with it? We may want to draw it, and for this we first need a 3D turtle, which is a geometric cursor, whose location and direction in space are always defined by its transformation matrix. To implement the recursive branching, we use a simple stack on which we'll push the current matrix whenever a new branch starts, and pop the previous one whenever a branch ends:

from numpy import *

class Turtle:

    def __init__(self):
        self.T = eye(4) # transformation matrix
        self.stack = []

    def getPos(self):
        return array([self.T[0][-1], self.T[1][-1], self.T[2][-1]])

    def getDir(self):
        return dot(self.T, [0, 1, 0, 0])[:3]

    def rotateX(self, a):
        self.T = dot(self.T, [[1, 0,      0,       0],
                              [0, cos(a), -sin(a), 0],
                              [0, sin(a), cos(a),  0],
                              [0, 0,      0,       1]])

    def rotateY(self, a):
        self.T = dot(self.T, [[cos(a),  0, sin(a), 0],
                              [0,       1, 0,      0],
                              [-sin(a), 0, cos(a), 0],
                              [0,       0, 0,      1]])

    def rotateZ(self, a):
        self.T = dot(self.T, [[cos(a), -sin(a), 0, 0],
                              [sin(a), cos(a),  0, 0],
                              [0,      0,       1, 0],
                              [0,      0,       0, 1]])

    def translate(self, d):
        self.T += [[0, 0, 0, d[0]],
                   [0, 0, 0, d[1]],
                   [0, 0, 0, d[2]],
                   [0, 0, 0, 0   ]]

    def forward(self, d=1):
        self.translate(self.getDir() * d)

    def pushMatrix(self):
        self.stack.append(self.T.copy()) # copy is important!

    def popMatrix(self):
        self.T = self.stack.pop()

The only missing piece is an "interpreter", and for this we can use the VPython (or Visual) 3D visualization module, whose primary quality is quite certainly.. ease of use. In fact the only thing we will ask this module is to draw a cylinder segment whenever we read an F character, as the rest can be handled by our turtle:

from visual import *

turtle = Turtle()
angle = radians(30)

for c in seq:

    if c == 'F':
        cylinder(pos=turtle.getPos(), axis=turtle.getDir(), radius=0.1)
        turtle.forward()

    elif c in '^v':
        turtle.rotateX(angle if c == 'v' else -angle)

    elif c in '+-':
        turtle.rotateY(angle if c == '+' else -angle)

    elif c in '<>':
        turtle.rotateZ(angle if c == '>' else -angle)

    elif c == '[':
        turtle.pushMatrix()

    elif c == ']':
        turtle.popMatrix()

With the addition of random angles and some parameters to control the evolution of segment length and radius, we can actually achieve interesting results like this:



Actually.. I cheated a little for this picture: I did not generate it using VPython, but rather with the Python bindings for VTK, a truly powerful toolkit, that I strongly recommend for any scientific visualization need.

There's a nice and freely available book about this very rich topic, co-written by its pioneer, Aristid Lindenmayer.