Hilbert curve indexing

Recently I wanted to generate hilbert curve indexing of 2D array as a part of one experiment. The experiment wasn't very successful itself, but I found the resulting Python code for curve generation prety neat to post it.

Here are a couple of interesting references on space-filling curve construction and theis applications:

  • http://blog.notdot.net/2009/11/Damn-Cool-Algorithms-Spatial-indexing-with-Quadtrees-and-Hilbert-Curves
  • http://www.win.tue.nl/~hermanh/stack/dagstuhl08-talk.pdf

So let's write some code!

In [2]:
def hilbert_curve(n):
    ''' Generate Hilbert curve indexing for (n, n) array. 'n' must be a power of two. '''
    # recursion base
    if n == 1:  
        return zeros((1, 1), int32)
    # make (n/2, n/2) index
    t = hilbert_curve(n//2)
    # flip it four times and add index offsets
    a = flipud(rot90(t))
    b = t + t.size
    c = t + t.size*2
    d = flipud(rot90(t, -1)) + t.size*3
    # and stack four tiles into resulting array
    return vstack(map(hstack, [[a, b], [d, c]]))   

That's it. Lets try some small n values.

In [3]:
print hilbert_curve(2)
[[0 1]
 [3 2]]
In [4]:
print hilbert_curve(4)
[[ 0  3  4  5]
 [ 1  2  7  6]
 [14 13  8  9]
 [15 12 11 10]]

Now let's show the generated indices with color. Note, how nearby cells tend to have close indices. That's why such curves are useful for spatial indexing.

In [5]:
_=imshow(hilbert_curve(64), interpolation='nearest')

Ok, time to plot the curve itself.

In [6]:
idx = hilbert_curve(32)
y, x = indices(idx.shape).reshape(2, -1)
x[idx.ravel()], y[idx.ravel()] = x.copy(), y.copy()

plot(x, y)
axis('equal')
axis('off')
_=ylim(ymin=-1)

This post in generated from IPython Notebook, which can be found in my GitHub repository

Comments !