Hilbert-Peano curve to scan image of arbitrary size

Question

I have written an implementation of Hilbert-Peano space filling curve in Python (from a Matlab one) to flatten my 2D image:

def hilbert_peano(n):
    if n<=0:
        x=0
        y=0
    else:
        [x0, y0] = hilbert_peano(n-1)
        x = (1/2) * np.array([-0.5+y0, -0.5+x0, 0.5+x0, 0.5-y0])
        y = (1/2) * np.array([-0.5+x0, 0.5+y0, 0.5+y0, -0.5-y0])

    return x,y

However, the classical Hilbert-Peano curve only works for multi-dimensionnal array whose shape is a power of two (ex: 256*256 or 512*512 in case of a 2D array (image)).

Does anybody know how to extend this to an array of arbitrary size?

jcerveny · Answer 1 · 2019-10-29T08:23:28.880

I had the same problem and have written an algorithm that generates a Hilbert-like curve for rectangles of arbitrary size in 2D and 3D. Example for 55x31: curve55x31

The idea is to recursively apply a Hilbert-like template but avoid odd sizes when halving the domain dimensions. If the dimensions happen to be powers of two, the classic Hilbert curve is generated.

def gilbert2d(x, y, ax, ay, bx, by):
    """
    Generalized Hilbert ('gilbert') space-filling curve for arbitrary-sized
    2D rectangular grids.
    """

    w = abs(ax + ay)
    h = abs(bx + by)

    (dax, day) = (sgn(ax), sgn(ay)) # unit major direction
    (dbx, dby) = (sgn(bx), sgn(by)) # unit orthogonal direction

    if h == 1:
        # trivial row fill
        for i in range(0, w):
            print x, y
            (x, y) = (x + dax, y + day)
        return

    if w == 1:
        # trivial column fill
        for i in range(0, h):
            print x, y
            (x, y) = (x + dbx, y + dby)
        return

    (ax2, ay2) = (ax/2, ay/2)
    (bx2, by2) = (bx/2, by/2)

    w2 = abs(ax2 + ay2)
    h2 = abs(bx2 + by2)

    if 2*w > 3*h:
        if (w2 % 2) and (w > 2):
            # prefer even steps
            (ax2, ay2) = (ax2 + dax, ay2 + day)

        # long case: split in two parts only
        gilbert2d(x, y, ax2, ay2, bx, by)
        gilbert2d(x+ax2, y+ay2, ax-ax2, ay-ay2, bx, by)

    else:
        if (h2 % 2) and (h > 2):
            # prefer even steps
            (bx2, by2) = (bx2 + dbx, by2 + dby)

        # standard case: one step up, one long horizontal, one step down
        gilbert2d(x, y, bx2, by2, ax2, ay2)
        gilbert2d(x+bx2, y+by2, ax, ay, bx-bx2, by-by2)
        gilbert2d(x+(ax-dax)+(bx2-dbx), y+(ay-day)+(by2-dby),
                 -bx2, -by2, -(ax-ax2), -(ay-ay2))

def main():
    width = int(sys.argv[1])
    height = int(sys.argv[2])

    if width >= height:
        gilbert2d(0, 0, width, 0, 0, height)
    else:
        gilbert2d(0, 0, 0, height, width, 0)

A 3D version and more documentation is available at https://github.com/jakubcerveny/gilbert

score 1 · Answer 2 · edited May 23 '17 at 12:14

1

I found this page by Lutz Tautenhahn:

"Draw A Space-Filling Curve of Arbitrary Size" (http://lutanho.net/pic2html/draw_sfc.html)

The algorithm doesn't have a name, he doesn't reference anyone else and the sketch suggests he came up with it himself.

I wonder if this is possible for a z order curve and how?

[1]Draw A Space-Filling Curve of Arbitrary Size

edited May 23 '17 at 12:14

Community

1
1

answered Jul 19 '16 at 16:08

Micromega

12,486
7
35
72

While its true the classical hilbert curve is based on power of 2 since it's a fractal you can also ignore at some cost. But to answer your question. No. – Micromega Jul 20 '16 at 07:44
How would you modify the original curve to get rid of the power of 2 hypothesis? – floflo29 Jul 20 '16 at 07:56
For example try a longer curve!? – Micromega Jul 20 '16 at 08:00
I thought about that but I will have to cut a part of the curve to fit my image and thus I will lose the property that two spatially close pixels will also be close on the space filling curve. – floflo29 Jul 20 '16 at 08:07

score 0 · Accepted Answer · answered Jul 21 '16 at 13:57

I finally choose, as suggested by Betterdev as adaptive curves are not that straigthforward [1], to compute a bigger curve and then get rid of coordinates which are outside my image shape:

# compute the needed order
order = np.max(np.ceil([np.log2(M), np.log2(N)]))
# Hilbert curve to scan a 2^order * 2^order image
x, y = hilbert_peano(order)
mat = np.zeros((2**order, 2**order))
# curve as a 2D array
mat[x, y] = np.arange(0, x.size, dtype=np.uint)
# clip the curve to the image shape
mat = mat[:M, :N]
# compute new indices (from 0 to M*N)
I = np.argsort(mat.flat)
x_new, y_new = np.meshgrid(np.arange(0, N, dtype=np.uint), np.arange(0, M, dtype=np.uint))
# apply the new order to the grid
x_new = x_new.flat[I]
y_new = y_new.flat[I]

[1] Zhang J., Kamata S. and Ueshige Y., "A Pseudo-Hilbert Scan Algorithm for Arbitrarily-Sized Rectangle Region"

Hilbert-Peano curve to scan image of arbitrary size

3 Answers3