Boids in QGIS part 2:
A note on testing

Posted on by Jim

An agent-based model typically consists of two components: an object describing our agent, and a control script that is responsible for running and tracking the system. If you want to think ecologically, this is a separation between our individual animals and their environment. Since we are ultimately going to use this agent in a QGIS plugin, the control script will be built into the QGIS plugin and I will cover that aspect later. For quick testing you can use this script. It creates a list of agents, and calls each of their updates and behaviour functions in a loop. The output is a csv file that you can import into QGIS to see the results as well as a timestamp that is readable by the QGIS TimeManager plugin (Graser et al.) This is part 2 of the series. Click here for part 1, or here for part 3.

Prerequisites:

We need to create a script in the same directory as our object and import both numpy and our object. I called my object boid.py and named the control script main.py

import numpy as np
from boid import Boid

Defining our study area:

Where do our Boids live? How far are they allowed to go? These numbers are in map units. Here I have chosen the Soquel Canyon MCA, a region at the centre of California’s Monterey Bay if you use UTM zone 15N as a CRS.

. . .
minX = 585076
minY = 4073066
maxX = 594870
maxY = 4078604

Creating our list of boids:

We are going to instantiate all of our boids in a Python list. We’ll start by putting them in a random location within the study area and decide on a number to instantiate. For now let’s go with 100.

. . .
# Figure out the size of our study area
xRange = maxX - minX
yRange = maxY - minY

# How many boids do we want?
nB = 100

# Let's generate the school!
school = [Boid(np.random.random() * 10 + minX + xRange/2, np.random.random() * 10 + minY + yRange/2, minX, minY, maxX, maxY) for _ in range(nB)]

Creating a file to store our results:

We want to output our results to a CSV file with our boid locations:

. . .
import numpy as np 
from boid import Boid

# Study area bounds (Soquel Canyon MPA in EPSG32610)
minX = 585076 
minY = 4073066 
maxX = 594870 
maxY = 4078604

# Figure out the size of our study area 
xRange = maxX - minX 
yRange = maxY - minY 

# How many boids do we want? 
nB = 100 

# Let's generate the school! 
school = [Boid(np.random.random() * xRange + minX, np.random.random() * yRange + minY, minX, minY, maxX, maxY) for _ in range(nB)]

# Create a file. Will give an error if it exists already 
f = open("positions.csv","xt") 

# Write our CSV header 
f.write("time,boid,x,y,direction\n")

# How many time steps do we want? 
timeSteps = 1000 

# Let's loop! 
for tS in range(timeSteps):
    for boid in school: 
        # Do the behaviours and update positions 
        boid.behave(school) 
        boid.update() 

        # Calculate the direction (optional but gives nice animations) 
        direction = np.arctan2(boid.velocity[0], boid.velocity[1]) * 57.29578 

        # Output our time stamp, id, position, and direction to the file 
        f.write(str(tS)) # Time Stamp 
        f.write(",") 
        f.write(str(boid.boidID)) # Unique boid ID 
        f.write(",") 
        f.write(str(boid.position[0])) # X 
        f.write(",") 
        f.write(str(boid.position[1])) # Y 
        f.write(",") 
        f.write(str(direction)) # Direction (see note) 
        f.write("\n") # Newline

Notes:

  • We’re creating a csv file called “Positions.csv” in the current directory. If Python can’t create the csv because it exists, it will error out, so you’ll need to delete or move this file between runs.
  • The direction is in degrees (as requested by QGIS) measured from the due east, with directions south of due east being expressed as negative numbers. Angles increase in the counterclockwise direction. You’ll need to define your symbology appropriately if you want to use it for symbol rotation. I like to use the GPS fish symbol, which normally points due west. To get this to play nicely, you can use direction + 90
  • The unique boid ID is generated in __init__ in our boid object. This is useful for tracking individuals, but if you use TimeManager, it may not be necessary to separate the symbols.
  • Time manager recognizes integer timestamps as Unix timestamps (seconds past the Unix Epoch) this means it will interpret all of our timestamps as being in the first 20 minutes of 1970. If this bothers you, you can add an arbitrary Unix timestamp to fix this behaviour, or be like me and turn off the timestamp display.

Wrapping it up:

This script is pretty close to the bare minimum necessary to run the simulation. The complete script looks like this:

import numpy as np 
from boid import Boid

# Study area bounds (Soquel Canyon MPA in EPSG32610)
minX = 585076 
minY = 4073066 
maxX = 594870 
maxY = 4078604

# Figure out the size of our study area 
xRange = maxX - minX 
yRange = maxY - minY 

# How many boids do we want? 
nB = 100 

# Let's generate the school! 
school = [Boid(*np.random.random() * xRange + minX, *np.random.random() * yRange + minY, minX, minY, maxX, maxY) for _ in range(nB)]

# Create a file. Will give an error if it exists already 
f = open("positions.csv","xt") 
# Write our CSV header 
f.write("time,boid,x,y,direction\n")

# How many time steps do we want? 
timeSteps = 1000 
# Let's loop! 
for tS in range(timeSteps):
  for boid in school: 
    # Do the behaviours and update positions 
    boid.behave(school) 
    boid.update() 
    # Calculate the direction (optional but gives nice animations) 
    direction = np.arctan2(boid.velocity[0], boid.velocity[1]) * 57.29578 
    # Output our time stamp, id, position, and direction to the file 
    f.write(str(tS)) # Time Stamp 
    f.write(",") 
    f.write(str(boid.boidID)) # Unique boid ID 
    f.write(",") 
    f.write(str(boid.position[0])) # X 
    f.write(",") 
    f.write(str(boid.position[1])) # Y 
    f.write(",") 
    f.write(str(direction)) # Direction (see note) 
    f.write("\n") # Newline