Thursday, November 20. 2014
This year's PostGIS day, I decided to celebrate with a little Conway's Game of Life fun inspired by Anita Graser's recent blog series Experiments with Game of Life.
The path I chose to simulate the Game of life is a little
different from Anita's. This variant exercises PostGIS 2.1+ raster mapalgebra and PostgreSQL 9.1+ recursive queries. Although you can do this with PostGIS 2.0, the map algebra syntax I am using is only supported in PostGIS 2.1+. My main disappointment is that because PostGIS does not yet support direct generation of animated gifs I had to settle for a comic strip I built unioning frames of rasters instead of motion picture. Hopefully some day my PostGIS animated gif dream will come true.
The Gosper's Glider Gun Initial state
Create table to store our initial state.
DROP TABLE IF EXISTS cellstates;
CREATE TABLE cellstates(id serial primary key, rast raster, state_name text);
We start with Gosper's Glider Gun initial state which is the smallest Cellular automaton Gun found. The gun can be represented as a matrix of 0s and 1s which with PostGIS ST_SetValues function (matrix variant), we can use to build our initial state raster (38x11 pixels).
INSERT INTO cellstates(rast, state_name)
WITH gunner AS ( SELECT
'{
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,1,0,0,0,0,0,0,0,0,0,0,0,0},
{0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0},
{0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,1,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0},
{0,1,1,0,0,0,0,0,0,0,0,1,0,0,0,0,0,1,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
{0,1,1,0,0,0,0,0,0,0,0,1,0,0,0,1,0,1,1,0,0,0,0,1,0,1,0,0,0,0,0,0,0,0,0,0,0,0},
{0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,1,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0},
{0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
{0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0}
}'::int[][] As initval )
SELECT ST_SetValues(
ST_AddBand(
ST_MakeEmptyRaster(38, 11, 0, 0, 1, -1,0,0,0),
1, '8BUI', 0, 0), 1, 1,1,
gunner.initval), 'gunner'::text As state_name
FROM gunner;
Note: that I'm using 8BUI instead of the more appropriate '1BB' pixel type. This is just so I can convert to PNG and other formats without hassle.
Game of Life for Map Algebra
In order to play this game, we'll need a map algebra callback function that simulates the hand of God in deciding about matters of life and death. I chose plpgsql to build the call back. I suspect the map algebra function would be much faster in PL/V8 but I'll save that for another day to experiment with.
Just as a refresher the rules of the game are:
- Any live cell with less than two live neighbours dies, as if caused by under-population.
- Any live cell with two or three live neighbours lives on to the next generation.
- Any live cell with more than three live neighbours dies, as if by overcrowding.
- Any dead cell with exactly three live neighbours becomes a live cell, as if by reproduction
Our call back function given the neighborhood matrix, computes the new state of the center pixel and returns that as the new value for that pixel.
CREATE OR REPLACE FUNCTION postgis_game_of_life_4MA(
value double precision[][][], pos integer[][],
VARIADIC userargs text[] DEFAULT NULL::text[])
RETURNS double precision
AS $$
DECLARE
nAlive int; z int; zRefCell int; yRefCell int; xRefCell int;
valRefCell double precision; var_result double precision;
BEGIN
nAlive := 0;
z := 1;
yRefCell := ( array_upper(value,2) - 1)/2 + 1;
xRefCell := ( array_upper(value,3) - 1)/2 + 1;
zRefCell := z;
valRefCell := COALESCE(value[zRefCell][yRefCell][xRefCell],0);
FOR y IN array_lower(value, 2)..array_upper(value, 2) LOOP
FOR x IN array_lower(value, 3)..array_upper(value, 3) LOOP
IF NOT ( (x,y,z) = (xrefCell, yrefCell,zRefcell)
OR COALESCE(value[z][y][x],0) = 0 ) THEN
nAlive = nAlive + 1;
END IF;
END LOOP;
END LOOP;
IF ( valRefCell = 1 AND (nAlive < 2 or nAlive > 3) ) THEN
var_result := 0;
ELSIF ( valRefCell = 0 AND nAlive = 3 ) THEN
var_result := 1;
ELSE
var_result := valRefCell;
END IF;
RETURN var_result;
END;
$$ LANGUAGE 'plpgsql' IMMUTABLE;
Let's play the Game of life
Now that we have all the pieces, we'll generate a reel of states by using the PostGIS ST_MapAlgebra function to drive our game of life call back function.
Note that the Game of Life is about 8 neighbors which translates in PostGIS neighborhood speak to x-distance=1 and y-distance=1 for a matrix that is 3x3.
The recursive call is a quickie way of repeating the same process on the newly created map algebraed raster. If you feel uncomfortable with recursive queries, or just want to black box this
you can replace it with a recursive plpgsql function. On my windows desktop this took about 200ms to generate the 10 frames.
DROP TABLE IF EXISTS gol_states;
WITH RECURSIVE
gol AS (
SELECT rast, 0 As nstate
FROM cellstates where state_name = 'gunner'),
gol2 AS (
SELECT rast, nstate
FROM gol
UNION ALL
SELECT ST_MapAlgebra (
gol2.rast
, 1,
'postgis_game_of_life_4MA(double precision[][][], integer[][],
text[])'::regprocedure,
'8BUI', 'FIRST', NULL, 1, 1
) As rast, gol2.nstate + 1 As nstate
FROM gol2
WHERE gol2.nstate < 10
)
SELECT nstate, rast INTO gol_states FROM gol2;
The above code creates a table gol_states where each row is a raster representing each state of the game from 0 to 10.
Since I couldn't create an animated gif, I opted to union all these rasters together so I have a single raster representing all the states.
The following code generates a single PNG with 20 pixel vertical/horizontal padding to separate each state and states running from left to right and down.
SELECT ST_AsPNG(
ST_Resize(
rast,ST_Width(rast)*4, ST_Height(rast)*4)
) As png
FROM (
SELECT ST_Union(
ST_SetUpperLeft(rast,
(nstate % 3)*(ST_Width(rast) + 20),
-(ST_Height(rast) + 20)*(nstate/3)::integer
) ) As rast
FROM gol_states) As reel;
To output the image you can use PSQL as described Raster Output PSQL
or you can skip the ST_ASPNG call, put out a dummy id and use QGIS Db Manager. I used my quickie NodeJS PostGIS viewer. The output of a 10 state run is:
|