Converting between Spiral Honeycomb Mosaic and axial hex coordinates

Question

There seems to be quite a lot of logic in dealing with hexagonal grids using a system of axial coordinates like the one described here.

In doing some reading about hex grids, though, I ran across the notion of the Spiral Honeycomb Mosaic (the links on that mirror are broken; the content is copied from here, about halfway down the page). It seems to be used mostly in image processing, but it's interesting from the standpoint that it assigns each hex a unique, sequential base-7 number as its identifier.

I'm not sure that it's actually useful for anything, but I'm wondering: is it possible to convert back and forth between the SHM numbering system and axial coordinates?

Both systems have methods to convert from a hex to the Cartesian coordinates of its center, so presumably one could calculate that out of one system, and then into the other. But that's messy and inefficient, and I'm convinced it's possible to do better.

Because of the fractal nature of the SHM, the conversion routine in the sample code provided relies on a loop; it seems like the general solution might involve a series?

Answer 1

Thanks for a fascinating puzzle! Yes, it looks like we can do better than a conversion through cartesian coordinates of hexagon centers. It can be done entirely with integer math, though I've included a rational in a matrix below to keep the notation concise.

You're right that both the encoding and decoding processes require loops. Fortunately, because SHM works by orders of magnitude, the loops will never need more than Log_7(n + 1) iterations, where n is the number of tiles in the grid, which is very reasonable.

To keep my notation symmetric, I've opted to use a 3-coordinate grid like the one below (what your example from Red Blob Games calls "Cube Coordinates"). If you prefer, you can omit the z coordinate and replace all instances of .z with -(x + y) to get "Axial Coordinates" instead.

For speed and cleanliness of code, I recommend using lookup tables for encoding and decoding. First, for decoding SHM into cube coordinates, you can build the table with these recurrences:

Rotate(v) = (-v.z, -v.x, -v.y)

Decode[,] = 2D array of ordered triples, [7 x (max_order_of_magnitude + 1)]

Decode[0, k] = (0, 0, 0) for all k = 0...max_order_of_magnitude
Decode[1, 0] = (0, -1, 1)

Decode[d + 1, k] = Rotate(Decode(d, k)) for d = 2...6, all k
Decode[d, k + 1] = Decode(d, k) + 2 * Rotate(Decode(d, k)) for all d, k

Then decoding becomes extremely simple:

SHMToPoint(code)
{
    point = (0, 0, 0)
    order = 0    

    while(code > 0)
    {
        digit = code % 7
        point += Decode(digit, order)
        code = floor(code/7)
        order++
    }
    return point
}

This works because the first row of the table stores the offsets of the cells encoded by the digits 0, 1, 2, 3, 4, 5, 6. The next row stores the offsets for 0, 10, 20, 30, 40, 50, 60, and so on. By summing the offsets contributed by each order of magnitude, you can decode any point within the range of the table.

For encoding, I think the following method will work, but I'll confess I haven't tested it. It uses a much simpler lookup table:

Encode = [0, 5, 1, 6, 3, 4, 2]

and a rotation & scaling matrix:

Untwist = | 5  -4   2 |
          | 2   5  -4 |
          |-4   2   5 | * 1/21

To encode a cube coordinate into SHM, then...

PointToSHM(point)
{
    code = 0
    magnitude = 1
    while(point != (0, 0, 0))
    {
        digitIndex = (((point.x - 2 * point.y) % 7) + 7) % 7
        code += magnitude *Encode[digitIndex]
        point = Untwist * (point - Decode[digit, 0])
        // treating "point" as a column vector, and rounding if necessary
        magnitude *= 7
    }
    return code
}

Note that the code generated here (and decoded above) is in standard numeric format (ie. the cell marked "10" in the SHM diagram will be assigned a code of 7 in base 10, or 111 in binary). If you want to display it in base 7, you'll need to use the appropriate radix when converting it to a string.

This works by figuring out the least significant digit of the SHM code at each iteration, (using the fact that these digits repeat in a way that tiles the plane regularly - the LUT is the sequence of digits read from (0, 0, 0) to (6, 0, -6), which repeats at an offset every row). Next it effectively zeroes that digit by subtracting the corresponding offset from our Decode table earlier, and scales the fractal down one level using the Untwist matrix. This acts like a right shift in base-7, in geometric form - it brings the next digit into the least significant position, bringing cell "10" to "1" and so on, where we can use the same trick again to read that digit on the next loop.

I hope that works for you!