Ch. 11 – Exploring Powerful Ideas

Topics:  Fractals, recursion, Fibonacci numbers, the Golden Ratio, Zipf’s Law, top-down design, Python dictionaries, defining Python classes, Python exceptions, animation, color gradients, Python complex numbers, cymatics and dynamical systems (boids).

In this chapter, we present several advanced computational concepts and algorithmic techniques. For some of them, we will show directly how to use them to make music (as we have done throughout the chapters so far). For others, we will leave it up to you to add a musical dimension to the sample code. More information is provided in the reference textbook.

Here is code from this chapter.

Fractals and Recursion

Fractals are self-similar objects (i.e., objects consisting of multiple parts, with smaller parts being reduced-size replicas of the larger parts).

In this section we explore a programming technique called, recursion.  Recursion mimics the subdivision process found in nature, e.g., the branching (starting with a single sprout) that results in a pine tree, etc.  Recursion is especially suited for creating fractals.

Fibonacci numbers and the Golden Ratio

Fibonacci numbers appear in many natural objects and processes, including seashells (see figure below), branching plants, and pine cones, among others. They also appear in the shape of hurricanes, and galaxies.

The code sample below (Ch. 11, p. 342) demonstrates one way to calculate the Fibonacci numbers through the process of recursion:

# Find the nth Fibonacci number using recursion.

def fib(n):

   if n == 0:    # simple case 1
      result = 0
   elif n == 1:  # simple case 2
      result = 1
   else:         # recursive case
      result = fib(n-2) + fib(n-1)

   return result

# now, let's test it
for i in range(10):
   print "fib(" + str(i) + ") =", fib(i)

Fractals and generating a golden tree

The golden ratio is found in natural and human-made artifacts. It is also found in the human body (e.g., the bones of our hands, the cochlea in our ears). It is also found in various musical works by Bach, Beethoven, Mozart, and others.

The following program (Ch. 11, p. 344) demonstrates how to generate a fractal via recursion.  This particular fractal, known as the Golden Tree, incorporates the golden ratio.

Here is the code:

# Demonstrates how to draw a golden tree using recursion.

from gui import *
from math import *

# create display
d = Display("Golden Tree", 250, 250)

# calculate phi to the highest accuracy Python allows
phi = (sqrt(5) - 1) / 2     # approx. 0.618033988749895

# recursive drawing parameters
depth = 10                  # amount of detail (or branching)
rotation = radians(60)      # branch angle is 60 degrees (need radians)
scale = phi                 # scaling factor of branches

# initial parameters
angle = radians(90)         # starting orientation is North
length = d.getHeight() / 3  # length of initial branch (trunk)
startX = d.getWidth() / 2   # start at bottom center
startY = d.getHeight() - 33

# recursive function for drawing tree
def drawTree(x, y, length, angle, depth):
   Recursively draws a tree of depth 'depth' starting at 'x', 'y'.

   global d, scale, rotation 

   #print "depth =", depth, "x =", x, " y =", y, " length =", length,
   #print "angle =", degrees(angle)

   # draw this line
   newX = x + length * cos( angle )    # calculate run
   newY = y - length * sin( angle )    # calculate rise
   d.drawLine(int(x), int(y), int(newX), int(newY))

   # check if we need more detail
   if depth > 1:

      # draw left branch - use line with length scaled by phi,
      # rotated counter-clockwise
      drawTree(newX, newY, length*phi, angle - rotation, depth-1)

      # draw right branch - use line with length scaled by phi,
      # rotated clockwise
      drawTree(newX, newY, length*phi, angle + rotation, depth-1)

# draw complete tree (recursively)
drawTree(startX, startY, length, angle, depth)

It generates this shape::

The Sierpinski triangle

Another fractal shape is the Sierpinski triangle. This fractal triangle consists of three smaller triangles (top, bottom left, bottom right) that have the same shape as the main one (see figure below).

Here is the code:

# Demonstrates how to draw the Sierpinski triangle using recursion.

from gui import *

# create display
d = Display("Sierpinski Triangle", 250, 250)

# recursive drawing parameters
depth = 6      # amount of detail (number of subdivisions) 

# initial points
x      = d.getWidth() / 2    # top corner
y      = 20
width  = d.getWidth() - 40
height = d.getHeight() - y - 30 

# recursive function for drawing triangle
def drawTriangle(x, y, width, height, depth):
   Recursively draws a Sierpinski triangle of depth 'depth' with top at 'x', 'y', and
   provided 'width' and 'height'.
   global d  # display to draw tree

   # get corners
   x1, y1 = x, y                     # top
   x2, y2 = x - width/2, y + height  # left
   x3, y3 = x + width/2, y + height  # right

   # reached the pixel level?
   if depth == 1:

      d.drawPolygon([x1, x2, x3], [y1, y2, y3])   # yes, so draw a triangle

   else:   # no, so continue subdividing

      # get top corners of subtriangles
      topX, topY     = x, y
      leftX, leftY   = x - width/4, y + height/2
      rightX, rightY = x + width/4, y + height/2

      # and draw them recursively
      drawTriangle(topX, topY, width/2, height/2, depth-1)      # draw top subtriangle
      drawTriangle(leftX, leftY, width/2, height/2, depth-1)    # draw left subtriangle
      drawTriangle(rightX, rightY, width/2, height/2, depth-1)  # draw right subtriangle

# draw complete Sierpinski triangle (recursively)
drawTriangle(x, y, width, height, depth)

It generates this shape::

Measuring Zipf proportions in a MIDI file

The following program (Ch. 11, p. 353) demonstrates how to read in a sequence of MIDI files and measure how well the note pitches in each song follow Zipf’s law. It is relatively easy to expand it to measure additional proportions, such as note durations.

It uses the following MIDI files:

Make sure they are saved in the same folder as the program below, prior to running it.

Here is the code:

# Demonstrates how to calculate Zipf metrics from MIDI files for
# comparative analysis. It calculates Zipf slopes and R^2 values
# for pitch.
# It also demonstrates how to use Python dictionaries to store
# collections of related items.
# Finally, it demonstrates how to implement an algorithm in a top-down
# fashion. First function encountered in the program performs the
# highest-level tasks, and any sub-tasks are relegated to lower-level
# functions.

from music import *
from zipf import *

# list of MIDI files to analyze
pieces = ["sonifyBiosignals.mid", "ArvoPart.CantusInMemoriam.mid",
          "Pierre Cage.Structures pour deux chances.mid"]

# define main function
def main( pieces ):
   """Calculates and outputs Zipf statistics of all 'pieces'."""

   # read MIDI files and count pitches
   for piece in pieces:

      # read this MIDI file into a score
      score = Score() # create an empty score
      Read.midi( score, piece ) # and read MIDI file into it

      # count the score's pitches
      histogram = countPitches( score ) 

      # calculate Zipf slope and R^2 value
      counts = histogram.values()
      slope, r2, yint = byRank(counts)

      # output results
      print "Zipf slope is", round(slope, 4), ", R^2 is", round(r2, 2),
      print "for", piece

   # now, all the MIDI files have been read into dictionary
   print # output one more newline

def countPitches( score ):
   """Returns count of how many times each individual pitch appears
   in 'score'.

   histogram = {} # holds each of the pitches found and its count

   # iterate through every part, and for every part through every
   # phrase, and for every phrase through every note (via nested
   # loops)
   for part in score.getPartList(): # for every part
      for phrase in part.getPhraseList(): # for every phrase in part
         for note in phrase.getNoteList(): # for every note in phrase

            pitch = note.getPitch() # get this note's pitch

            # count this pitch, if not a rest
            if (pitch != REST):

               # increment this pitch's count (or initialize to 1)
               histogram[ pitch ] = histogram.get(pitch, 0) + 1

         # now, all the notes in this phrase have been counted
      # now, all the phrases in this part have been counted
   # now, all the parts have been counted

   # done, so return counts
   return histogram

# start the program
main( pieces )

Implementing the Note class

The following program (Ch. 11, p. 361demonstrates how to define a simplified version of the Note class. This class is used to store related information about musical notes (i.e., pitch, duration, dynamic, and panning). It also provides functions to retrieve (get) and update (set) the information encapsulated in a Note object.

Here is the code:

# It demonstrates how to create a class to encapsulate a musical
# note.  (This is a simplified version of Note class found in the
# music library.)

from music import *

class Note:

   def __init__(self, pitch=C4, duration=QN, dynamic=85, panning=0.5):
      """Initializes a Note object."""

      self.pitch = pitch        # holds note pitch (0-127)
      self.duration = duration  # holds note duration (QN = 1.0)
      self.dynamic = dynamic    # holds note dynamic (0-127)
      self.panning = panning    # holds note panning (0.0-1.0)

   def getPitch(self):
      """Returns the note's pitch (0-127)."""

      return self.pitch

   def setPitch(self, pitch):
      """Sets the note's pitch (0-127)."""

      # first ensure data integrity, then update
      if 0 <= pitch <= 127:  # is pitch in the right range?
         self.pitch = pitch     # yes, so update value
      else:                  # otherwise let them know
         print "TypeError: Note.setPitch(): pitch ranges from",
         print "0 to 127 (got " + str(pitch) + ")"
         #raise TypeError("Note.setPitch(): pitch ranges from 0 to 127 (got " + str(pitch) + ")")

   def getDuration(self):
      """Returns the note's duration."""

      return self.duration

   def setDuration(self, duration):
      """Sets the note's duration (a float)."""

      if type(duration) == type(1.0):    # is duration a float?
         self.duration = float(duration)    # yes, so update value
      else:                  # otherwise let them know
         print "TypeError: Note.setDuration(): duration must be",
         print "a float (got " + str(duration) + ")"
         # raise TypeError("Note.setDuration(): duration must be a float (got " + str(duration) + ")")

if __name__ == '__main__':

   # create a note
   print "n = Note(C4, QN)"
   n = Note(C4, QN)
   print "n.getPitch() =", n.getPitch()
   print "n.getDuration() =", n.getDuration()

Implementing a Slider control

The following program (Ch. 11, p. 369) demonstrates how to define a SliderControl class, which can be used to create individual sliders.

A SliderControl object has a display with a single slider. This can be used to control the value of any program variable interactively (as demonstrated below).

Here is the code:

# It creates a simple slider control surface.  It allows arbitrary
# numeric values for max, min, and start, as long as
# min < start < max.  However, the accuracy (or resolution) of the
# slider control (its ability to discriminate between two
# "consecutive" values) depends on screen resolution and mouse
# accuracy. (If the min to max range has too many values in it
# (e.g., 1 to 1000), some will not be "reachable".)
from gui import *

class SliderControl:

   def __init__(self, title="Slider", updateFunction=None,
                      minValue=10, maxValue=1000, startValue=None,
                      x=0, y=0):
      """Initializes a SliderControl object."""

      # holds the title of the control display
      self.title = title    

      # external function to call when slider is updated
      self.updateFunction = updateFunction  

      # determine slider start value
      if startValue == None:  # if startValue is undefined

         # start at middle point between min and max value
         startValue = (minValue + maxValue) / 2   

      # create slider
      self.slider = Slider(HORIZONTAL, minValue, maxValue, startValue, self.setValue)

      # create control surface display
      self.display = Display(self.title, 250, 50, x, y)

      # add slider to display
      self.display.add(self.slider, 25, 10)

      # finally, initialize value and title (using 'startValue')
      self.setValue( startValue ) 

   def setValue(self, value):
      """Updates the display title, and calls the external
         update function with given 'value'.

      self.display.setTitle(self.title + " = " + str(value))

if __name__ == '__main__':

   # one example
   def printValue(value):
      print value

   control = SliderControl("Print", printValue)  # just the title and update function, all else using default values

   # another example (demonstrating how to change a global variable - more useful)
   testValue = -1

   def changeAndPrintValue(value):
      global testValue   

      testValue = value
      print testValue

   control1 = SliderControl("Change and Print", changeAndPrintValue, 0, 100, 50, 100, 100)  # modify program variable

   # and another one
   from music import *    # for Play.note()

   def playNote(pitch):
      """Plays an 1-sec note with the provided pitch."""
      Play.note(pitch, 0, 1000)

   control2 = SliderControl("Play Note", playNote, 0, 100, 50, 200, 200)  # modify program variable

When run by itself, it generates the following sliders:

Animation – a revolving music sphere

This code sample (Ch. 11, p. 375) demonstrates how to create an animation involving a revolving sphere.

We model the sphere by drawing points on a GUI display. As the sphere rotates, it plays musical notes. Rotation speed is adjusted using a slider on a separate display.

Here is the code:

# Demonstrates how to create an animation of a 3D sphere using regular points
# on a Display.  The sphere is modeled using points on a spherical coordinate
# system (see
# It converts from spherical 3D coordinates to cartesian 2D coordinates (mapping
# the z axis to color).  When a point passes the primary meridian (the imaginary
# vertical line closest to the viewer), a note is played based on its latitude
# (low to high pitch).  Also the point turns red momentarily.
# Based on code by Uri Wilensky (1998), distributed with NetLogo
# under a Creative Commons Attribution-NonCommercial-ShareAlike 3.0 License.
from gui import *
from music import *
from random import *
from math import *
from sliderControl import *  

class MusicalSphere:
   """Creates a revolving sphere that plays music."""

   def __init__(self, radius, density, velocity=0.01, frameRate=30):
      Construct a revolving sphere with given 'radius', 'density'
      number of points (all on the surface), moving with 'velocity' angular
      (theta / azimuthal) velocity, at 'frameRate' frames (or movements) per
      second.  Each point plays a note when crossing the zero meridian (the
      sphere's meridian (vertical line) closest to the viewer).

      # musical parameters
      self.instrument = PIANO
      self.scale = MAJOR_SCALE
      self.lowPitch = C1
      self.highPitch = C6
      self.noteDuration = 2000    # milliseconds (2 seconds)  

      Play.setInstrument(self.instrument, 0)   # set the instrument

      # visual parameters
      self.display = Display("3D Sphere", radius*3, radius*3)  # display to draw sphere
      self.display.setColor( Color.BLACK )  # set background color to black

      self.radius = radius       # how wide circle is
      self.numPoints = density   # how many points to draw on sphere surface
      self.velocity = velocity   # how far it rotates per animation frame
      self.frameRate = frameRate # how many animation frames to do per second

      self.xCenter = self.display.getWidth() / 2   # to place circle at display's center
      self.yCenter = self.display.getHeight() / 2

      # sphere data structure (parallel lists)
      self.points      = []    # holds the points
      self.thetaValues = []    # holds the points' rotation (azimuthal angle)
      self.phiValues   = []    # holds the points' latitude (polar angle)

      # timer to drive animation
      delay = 1000 / frameRate   # convert from frame rate to timer delay (in milliseconds)
      self.timer = Timer(delay, self.movePoints)   # timer to schedule movement

      # control surface for animation frame rate
      xPosition = self.display.getWidth() / 3    # set initial position of display
      yPosition = self.display.getHeight() + 45
      self.control = SliderControl(title="Frame Rate", updateFunction=self.setFrameRate,
                                   minValue=1, maxValue=120, startValue=self.frameRate,
                                   x=xPosition, y=yPosition)

      # orange color gradient (used to display depth, the further away, the darker)
      black = [0, 0, 0]         # RGB values for black
      orange = [251, 147, 14]   # RGB values for orange
      white = [255, 255, 255]   # RGB values for white

      # create list of gradient colors from black to orange, and from orange to white
      # (a total of 25 colors)
      self.gradientColors = colorGradient(black, orange, 12) + colorGradient(orange, white, 12) + [white]  # remember to include the final color

      self.initSphere()      # create the circle

      self.start()           # and start rotating!

   def start(self):
      """Starts sphere animation."""

   def stop(self):
      """Stops sphere animation."""

   def setFrameRate(self, frameRate=30):
      """Controls speed of sphere animation (by setting how many times per second to move points)."""

      delay = 1000 / frameRate   # convert from frame rate to delay between each update (in milliseconds)

   def initSphere(self):
      """Generate a sphere of 'radius' out of points (placed on the surface of the sphere)."""

      for i in range(self.numPoints):     # create all the points

         # get random spherical coordinates for this point
         r = self.radius                                  # all points are placed *on* the surface
         theta = mapValue( random(), 0.0, 1.0, 0.0, 2*pi) # random rotation (azimuthal angle)
         phi = mapValue( random(), 0.0, 1.0, 0.0, pi)     # random latitude (polar angle)

         # project from spherical to cartesian 3D coordinates (z is depth)
         x, y, z = self.sphericalToCartesian(r, phi, theta)  

         # convert depth (z) to color
         color = self.depthToColor(z, self.radius)      

         # create a point (with thickness 1) at these x, y coordinates
         point = Point(x, y, color, 1)      

         # remember this point and its spherical coordinates (r equals self.radius for all points)
         # (append data for this point to the three parallel lists)
         self.points.append( point )
         self.phiValues.append( phi )
         self.thetaValues.append( theta )

         # add this point to the display
         self.display.add( point )

   def sphericalToCartesian(self, r, phi, theta):
      """Convert spherical to cartesian coordinates."""   

      # adjust rotation so that theta is 0 at max z (i.e., closest to viewer)
      x = int( r * sin(phi) * cos(theta + pi/2) )   # horizontal axis (pixels are int)
      y = int( r * cos(phi) )                       # vertical axis
      z = int( r * sin(phi) * sin(theta + pi/2) )   # depth axis      

      # move sphere's center to display's center
      x = x + self.xCenter
      y = y + self.yCenter

      return x, y, z

   def depthToColor(self, depth, radius):
      """Map 'depth' to color using the 'gradientColors' RGB colors."""

      # create color based on position (points further away have less luminosity)
      colorIndex = mapValue(depth, -self.radius, self.radius, 0, len(self.gradientColors))  # map depth to color index
      colorRGB = self.gradientColors[colorIndex]                    # get corresponding RBG value
      color = Color(colorRGB[0], colorRGB[1], colorRGB[2])          # and create the color

      return color

   def movePoints(self):
      """Rotate points on y axis as specified by angular velocity."""

      for i in range(self.numPoints):
         point = self.points[i]                                   # get this point
         theta = (self.thetaValues[i] + self.velocity) % (2*pi)   # increment angle to simulate rotation
         phi = self.phiValues[i]                                  # get latitude (altitude)

         # convert from spherical to cartesian 3D coordinates
         x, y, z = self.sphericalToCartesian(self.radius, phi, theta)

         if self.thetaValues[i] > theta:   # did we just cross the primary meridian?
            color = Color.RED                  # yes, so sparkle for a split second
            pitch = mapScale(phi, 0, pi, self.lowPitch, self.highPitch, self.scale)  # phi is latitude
            dynamic = randint(0, 127)          # random dynamic
            Play.note(pitch, 0, self.noteDuration, dynamic)  # and play a note (based on latitude)

         else:   # we are not on the primary meridian, so
            # convert depth to color (points further away have less luminosity)
            color = self.depthToColor(z, self.radius)

         # adjust this point's position and color
         self.display.move(point, x, y)

         # now, remember this point's new theta coordinate
         self.thetaValues[i] = theta

# create a sphere
sphere = MusicalSphere(radius=200, density=200, velocity=0.01, frameRate=30)

This code generates the following animation. Since the points are placed randomly, the music generated will be different every time:

Boids – implementing flocking behavior

The following program (Ch. 11, p. 389) demonstrates the types of behaviors that emerge from self-adapting agents (boids), which sense their surroundings and modify their behavior (movement) accordingly.  These boids (points) move around the display, staying in close proximity to other boids, yet avoiding them, while collectively moving towards the mouse (when the mouse is dragged).

This code may be adapted to resemble the movement of birds, fish, and other animal collectives.  There are numerous ways to add sound to it (proximity, cohesion, movement, speed of movement, rean collisions, and so on).  The possibilities are endless (so, none have been provided here). Explore.

Here is the code:

# This program simulates 2D boid behavior.
# See and

from gui import *
from math import *
from random import *

# universe parameters
universeWidth  = 1000   # how wide the display
universeHeight = 800    # how high the display

# boid generation parameters
numBoids   = 200        # from 2 to as much as your CPU can handle
boidRadius = 2          # radius of boids
boidColor  = Color.BLUE # color of boids

# boid distance parameters
minSeparation  = 30     # min comfortable distance between two boids
flockThershold = 100    # boids closer than this are in a local flock

# boid behavior parameters (higher means quicker/stronger)
separationFactor = 0.01 # how quickly to separate
alignmentFactor = 0.16  # how quickly to align with local flockmates
cohesionFactor = 0.01   # how quickly to converge to attraction point
frictionFactor = 1.1    # how hard it is to move (dampening factor)

### define boid universe, a place where boids exist and interact ####
class BoidUniverse:
   """This is the boid universe, where boids exist and interact.
      It is basically a GUI Display, with boids (moving Circles)
      added to it.  While boids are represented as circles, they
      have a little more logic to them - they can sense where other
      boids are, and act accordingly.  While individual boids have
      simple rules for sensing their environment and reacting,
      very intricate, complex, naturally-looking patterns of behavior
      emerges, similar to those of birds flying high in the sky
      (among others).  The rules of behavior are the same for all
      boids, and are defined in the Boid class (a sister class to
      this one).

   def __init__(self, title = "", width = 600, height = 400,

      self.display = Display(title, width, height) # universe display

      self.boids = []                              # list of boids

      # holds attraction point for boids (initially, universe center)
      self.attractPoint = complex(width/2, height/2)  

      # create timer
      delay = 1000 / frameRate     # convert frame rate to delay (ms)
      self.timer = Timer(int(delay), self.animate)       # animation timer

      # when mouse is dragged, call this function to set the
      # attraction point for boids
      self.display.onMouseDrag( self.moveAttractionPoint )

   def start(self):
      """Starts animation."""
      self.timer.start()   # start movement!

   def stop(self):
      """Stops animation."""
      self.timer.stop()    # stop movement

   def add(self, boid):
      """Adds another boid to the system."""

      self.boids.append( boid )          # remember this boid
      self.display.add( )    # add a circle for this boid

   def animate(self):
      """Makes boids come alive."""

      ### sensing and acting loop for all boids in the universe !!!
      for boid in self.boids:   # for every boid

         # first observe other boids and decide how to adjust movement
         boid.sense(self.boids, self.attractPoint)   

         # and then, make it so (move)!

   def moveAttractionPoint(self, x, y):
      """Update the attraction point for all boids."""
      self.attractPoint = complex(x, y)

### define the boids, individual agents who can sense and act #######
class Boid:
   """This a boid.  A boid is a simplified bird (or other species)
      that lives in a flock.  A boid is represented as a circle,
      however, it has a little more logic to it - it can sense where
      other boids are, and act, simply by adjusting its direction of
      movement.  The new direction is a combination of its reactions
      from individual rules of thumb (e.g., move towards the center
      of the universe, avoid collisions with other boids, fly in the
      same general direction as boids around you (follow the local
      flock, as you perceive it), and so on.  Out of these simple
      rules intricate, complex behavior emerges, similar to that of
      real birds (or other species) in nature.

   def __init__(self, x, y, radius, color,
                initVelocityX=1, initVelocityY=1 ):
      """Initialize boid's position, size, and initial
         velocity (x, y).

      # a boid is a filled circle = Circle(x, y, radius, color, True) 

      # set boid size, position
      self.radius = radius                # boid radius
      self.coordinates = complex(x, y)    # boid coordinates (x, y)

      # NOTE: We treat velocity in a simple way, i.e., as the
      # x, y displacement to add to the current boid coordinates,
      # to find where to move its circle next.  This moving is done
      # once per animation frame.

      # initialize boid velocity (x, y)
      self.velocity = complex(initVelocityX, initVelocityY)  

   def sense(self, boids, center):
      Sense other boids' positions, etc., and adjust velocity
      (i.e., the displacement of where to move next).

      # use individual rules of thumb, to decide where to move next

      # 1. Rule of Separation - move away from other flockmates
      #                         to avoid crowding them
      self.separation = self.rule1_Separation(boids)

      # 2. Rule of Alignment - move towards the average heading
      #                        of other flockmates
      self.alignment = self.rule2_Alignment(boids)

      # 3. Rule of Cohesion - move toward the center of the universe
      self.cohesion = self.rule3_Cohesion(boids, center)

      # 4. Rule of Avoidance: move to avoid any obstacles
      #self.avoidance = self.rule4_Avoidance(boids)

      # create composite behavior
      self.velocity = (self.velocity / frictionFactor) + \
                      self.separation + self.alignment + \

   def act(self, display):
      """Move boid to a new position using current velocity."""

      # Again, we treat velocity in a simple way, i.e., as the
      # x, y displacement to add to the current boid coordinates,
      # to find where to move its circle next.

      # update coordinates
      self.coordinates = self.coordinates + self.velocity

      # get boid (x, y) coordinates
      x = self.coordinates.real  # get the x part
      y = self.coordinates.imag  # get the y part

      # act (i.e., move boid to new position)
      display.move(, int(x), int(y) )

   ##### steering behaviors ####################
   def rule1_Separation(self, boids):
      """Return proper velocity to keep separate from other boids,
         i.e., avoid collisions.

      newVelocity = complex(0, 0)  # holds new velocity

      # get distance from every other boid in the flock, and as long
      # as we are too close for comfort, calculate direction to
      # move away (remember, velocity is just an x, y distance
      # to travel in the next animation/movement frame)
      for boid in boids:           # for each boid

         separation = self.distance(boid)   # how far are we?

         # too close for comfort (excluding ourself)?
         if separation < minSeparation and boid != self:
            # yes, so let's move away from this boid
            newVelocity = newVelocity - \
                          (boid.coordinates - self.coordinates)

      return newVelocity * separationFactor  # return new velocity

   def rule2_Alignment(self, boids):
      """Return proper velocity to move in the same general direction
         as local flockmates.
      totalVelocity = complex(0, 0) # holds sum of boid velocities
      numLocalFlockmates = 0        # holds count of local flockmates
      # iterate through all the boids looking for local flockmates,
      # and accumuate all their velocities
      for boid in boids:
         separation = self.distance(boid)    # get boid distance
         # if this a local flockmate, record its velocity
         if separation < flockThershold and boid != self:                      
            totalVelocity = totalVelocity + boid.velocity             
            numLocalFlockmates = numLocalFlockmates + 1       
      # average flock velocity (excluding ourselves)       
      if numLocalFlockmates > 0:
         avgVelocity = totalVelocity / numLocalFlockmates
         avgVelocity = totalVelocity
      # adjust velocity by how quickly we want to align
      newVelocity = avgVelocity - self.velocity
      return newVelocity * alignmentFactor  # return new velocity

   def rule3_Cohesion(self, boids, center):
      """Return proper velocity to bring us closer to center of the

      newVelocity = center - self.coordinates

      return newVelocity * cohesionFactor  # return new velocity

   ##### helper function ####################
   def distance(self, other):
      """Calculate the Euclidean distance between this and
         another boid.

      xDistance = (self.coordinates.real - other.coordinates.real)
      yDistance = (self.coordinates.imag - other.coordinates.imag)

      return sqrt( xDistance*xDistance + yDistance*yDistance )

# start boid simulation
universe = BoidUniverse(title="Boid Flocking Behavior",
                        width=universeWidth, height=universeHeight,

# create and place boids
for i in range(0, numBoids):

   # get random position for this boid
   x = randint(0, universeWidth)
   y = randint(0, universeHeight)

   # create a boid with random position and velocity
   boid = Boid(x, y, boidRadius, boidColor, 1, 1)
   universe.add( boid )

# animate boids

It generates an interactive animation similar to this:

In closing…

This website provides source code and related materials (e.g., API, etc.).

For more details and thorough explanations see the book. Also, the book offers observations, facts, and exercises for further musical explorations, and ways to create music via programming. Enjoy.

Ars longa, vita brevis…