Poisson-Disc Sampling and Generative Art

A little while ago I made a post about recreating some generative art I'd seen on the web by Espen Kluge and cam to a point in the code where I had to generate random points on the image. I didn't really think much of it at the time. However, this turns out to be a really interesting topic within the areas of game dev or generative art. How do you scatter points within an area that will be different every time but more evenly distributed across the plane? What I'd found out is that using the random function will not give you a truly random sampling. Some points will cluster together on the area, not a nice distribution across the image area.

Random That's not Truly Random

The short answer for this, is that the randomness used by Processing or P5 or Javascript is not actually a random process. It makes use of what's called a Psuedo-Random Number Generator. The distinction (which I learned about here and here. Essentially, the computer will use some internalized seed value to generate a number and the seed will change on every subsequent running of the random function. This means that if we knew the state of the random seed, then the random function would actually be predictable and determined.

In fact, processing has a function allowing you to set a seed value originally, such that the random function will give you the same output every time you run the function (refer here).

So random will actually give you a patterned output rather than a smooth distribution. This is where Poisson-Disc Sampling comes in. The technique in the algorithm is to split the area into a grid, keep track of which points you've laid down, and do so in O(n) time where n is the size of points you have. Pretty sick!

The Algorithm

I'm partly going to distill what Dan Shiffman goes over in his coding train video here and give you just the basics of the algoirthm.

The high level view of the algorithm is to section up the space into a grid, and each iteration will randomly choose a point closest to the last and check that this point is not within a certain distance to any other point on the grid. I'm gonna specifically do this in 2 dimensions but this can be extended into any arbitrary amount.

The Variables

width & height : How big the sampling area is. These are given to us for free in p5 and processing.

r : The minimum distance separating each point. Because we're evenly distributing our sampling, the algorithm must know how far apart each sampled point must be.

k : The max number of guesses the algorithm can make to place a point before moving on. This stops the algorithm from trying to place a point that is too close to nearby points.

grid : This is a one dimensional array that holds all the points of the space you sample on. Using nested for loops you'll be able to access the items in the array according to their position in the space (more on this below).

active: This is also a one dimensional array that holds all the points in the sample that have been placed. This will be handy for continuous point generation.

To the Code!

I'm gonna use processing and Java for this so the types of each variable will be:

import java.util.ArrayList;
float k = 30;
float r = 10;
PVector[] grid;
ArrayList<PVector> active = new ArrayList<PVector>();

The grid won't change significantly upon running it so there's no need to use the ArrayList data structure. However, the active list requires pushing and popping off the array so this needs to change throughout.

Step 1: Generate a random point in the grid

The algorithm kicks off by randomly locating a point in the sampling space and adding it to the active list. My code looked like this:

import java.util.ArrayList;
float k = 30;
float r = 10;
int cols;
int rows;
float w = r / sqrt(2);
PVector[] grid;
ArrayList<PVector> active = new ArrayList<PVector>();

void setup() {
  size(400,400);
  background(0);
  cols = floor(width / w);
  rows = floor(height / w);

  grid = new PVector[rows*cols];
  for (int i = 0; i < cols * rows; i++) {
    grid[i] = null;
  }

  PVector point = new PVector(random(width), random(height));

  int i = floor(point.x/w);
  int j = floor(point.y/w);

  grid[i + j * cols] = point;
  active.add(point);
}

Other than the normal processing setup stuff, I've initialized the amount of columns and rows that we need, created the amount of space we need in the grid by multiplying the cols by the rows and a variable w that will be the length of a circle of radius r that encapsulates a square. See:

This makes it impossible for two sampled points to be in the same grid cell. We initialize a random point vector in the space using processing and translate that point into a position on the grid and add that point into our active points list.

Step 2: Attempt to Place a New Point

Now is the trickiest part of the algorithm. We loop take a sample from the active array and try to generate a new point that is at least r distance away but less than 2 * r. We'll do this k amount of times so that we aren't stuck in an infinite loop. Here's the code I wrote to accomplish this:

void draw() {
  background(0);

  if (active.size() > 0) {
    int i = floor(random(active.size()));
    PVector pos = active.get(i);
    for (int j = 0; j < k; j++) {
      PVector sample = PVector.random2D();
      float m = random(r, 2 * r);
      sample.setMag(m);
      sample.add(pos);
      if (testSample(sample) == true) {
        active.add(sample);
        int x = floor(sample.x / w);
        int y = floor(sample.y / w);
        grid[x + y * cols] = sample;
        break;
      } else if (j == k - 1) {
        active.remove(i);
      }
    }
  }
}

Boolean testSample(PVector sample) {
  int col = floor(sample.x / w);
  int row = floor(sample.y / w);
  //println(col, row, cols, rows, grid[col + row * cols]);
  if (col > 0 && row > 0 && col < cols - 1 && row < rows - 1 && grid[col + row * cols] == null) {
    for (int i = -1; i <= 1; i++) {
      for (int j = -1; j <= 1; j++) {
        int index = (col + i) + (row + j) * cols;
        PVector neighbor = grid[index];
        if (neighbor != null) {
          float d = PVector.dist(sample, neighbor);
          if (d < r) {
            return false;
          }
        }
      }
    }
    return true;
  }
  return false;
}

I'll start from the top and move my way down now. So since the draw loop is run over and over, we can use that as a while loop. So if the active array is empty, we have no position to generate samples from, which means we'd have generated everything. Next, we'll randomly grab an element in the active array. We'll randomly make a 2D vector, set it's magnitude or length to between r and 2*r, and then add the element we're generating around to this new vector. This is in part due to nice vector attributes.

Once we've gotten our generated vector, we have to test whether this vector is not within r distance to another point. This takes us to the method "testSample" that I've written. It takes the sample we've made and checks all the adjacent grid locations around it to see if it's too close to them. It's important to note that not having the grid would mean we'd have to check all the points we've generated so far. Assigning them to grid locations means we can quickly check.

We'll loop through going between one column above and below our sample as well as one row above and below our sample. There was some weird error handling that had to be done if we were at the first and last row and column, and if that grid point had already been generated. Finally, if we encounter a neighbor and that neighbor is too close to our sample, the method returns false. And if we haver checked all adjacent grid cells and no red flags, then this point is good and we can return true.

if (testSample(sample) == true) {
        active.add(sample);
        int x = floor(sample.x / w);
        int y = floor(sample.y / w);
        grid[x + y * cols] = sample;
        break;
      } else if (j == k - 1) {
        active.remove(i);
      }

So if this tested sample is good, we find it's location in the grid, add the point to our grid and add it into the active list because it can then generate a further. However, if we couldn't generate a point AND the the loop variable j is about to break out of the loop (i.e. we've generated k points) then we pop the point we've been using to generate samples because it couldn't in k trials.

And BOOM, we've gotten a full poisson-disc sampling algorithm ready to go. If you wanted to take this out of processing, just replace the "if (active.size() > 1)" with a while loop and it should work just fine.