Morton Key Function for PostgreSQL/PostGIS

In my article on the creation and use of a space keys for the spatial sorting of spatial data, I show how to implement a Morton key for PostgreSQL/PostGIS pl/pgSQL.

The function is as follows:

 CREATE OR REPLACE FUNCTION ST_Morton (p_col int8, p_row int8)
  RETURNS int8
 AS
 $$
 /*  this procedure calculates the Morton number of a cell
     at the given row and col[umn]
     Written:  D.M. Mark, Jan 1984;
     Converted to Vax/VMS: July 1985
     Converted to PostgreSQL, Simon Greener, 2010
 */
 DECLARE
    v_row          int8 := 0;
    v_col          int8 := 0;
    v_key          int8;
    v_level        int8;
    v_left_bit     int8;
    v_right_bit    int8;
    v_quadrant     int8;
 BEGIN
    v_row   := p_row;
    v_col   := p_col;
    v_key   := 0;
    v_level := 0;
    WHILE ((v_row>0) OR (v_col>0)) LOOP
      /* Split off the row (left_bit) and column (right_bit) bits and
      then combine them to form a bit-pair representing the quadrant */
      v_left_bit  := v_row % 2;
      v_right_bit := v_col % 2;
      v_quadrant  := v_right_bit + 2*v_left_bit;
      v_key       := v_key + ( v_quadrant << (2*v_level) );
      /* row, column, and level are then modified before the loop continues */
      v_row := v_row / 2;
      v_col := v_col / 2;
      v_level := v_level + 1;
    END LOOP;
    RETURN (v_key);
 END;
 $$
   LANGUAGE 'plpgsql' IMMUTABLE;

Now, let’s test it with some simple SQL:

 SELECT *
   FROM (SELECT a.gcol,
                b.grow,
                st_morton( a.gcol, b.grow ) AS MortonKey,
                ST_MakeBox2D(ST_MakePoint(a.gcol, b.grow)::geometry,
                             ST_MakePoint(a.gcol+100, b.grow+100)::geometry)::geometry AS geometry           FROM (SELECT 0 + g AS gcol FROM generate_series(0,7,1) AS g) AS a,
                (SELECT 0 + g AS grow FROM generate_series(0,7,1) AS g) AS b
         ) AS foo
  ORDER BY mortonkey, gcol;

Results

I hope this is of use to someone.