Step 1: Don’t store depths, store the quadtree

The depths array is a pre-order listing of a quadtree, and a quadtree can be written down with a single yes/no question per node: does this cell split into four? For the Amsterdam and Diemen example from earlier, walking the tree from the root:

That is 13 bits for what used to be 128 bits (16 bytes) of depths. A cell at the maximum depth cannot split, so it does not even get a question. Even for this small example, that is already a reduction of 9.8×.

Those 13 bits are still not stored as 13 bits. Each question is asked with some knowledge of what was already decoded: the depth of the node, and whether the records to the west and south of it are deeper than it. If your neighbours had to split to follow a border, there is a good chance you have to as well. If your neighbours are big cells, you are probably in the middle of a city and you are probably a big cell too.

Step 2: Guess the palette index from the neighbours

This is the part I did not see coming. In my head, I was thinking of run length encoding, which is just adding counters after every data byte to mark how often they happen. But in Morton order, the record decoded just before you is often not even touching you: the Z-curve jumps around. What is always decoded before you are the cells directly to your west and south, because they have a lower Morton code.

So every leaf asks, in order:

1. Is it the same place as my west neighbour? In the world file, this is yes for 74% of the records.
2. If not: is it the same place as my south neighbour? Yes for another 19%.
3. If not: which one is it? (the remainder 7% of the cases) This is coded as a position in a move-to-front list of palette indices, where the most recently used index is at the front. The west and south answers are skipped in that list, since we already know it is neither of them. The position is written in Elias gamma, so small numbers are short.

The decoder finds the neighbours cheaply: for every row and every column of the chunk, it remembers the last record it decoded there. In Morton order, the last record in your row is the one directly to your west.

Step 3: A range coder instead of zlib

After steps 1 and 2, a chunk is nothing but a long list of yes/no decisions. Those bits are then fed into a range coder, the same kind as the one inside LZMA (and so in xz and 7-Zip). A range coder is a form of arithmetic coding: every decision comes with a probability, and a decision that was very likely costs very little, less than a bit:

The probabilities are not stored anywhere. Every one starts at 50% and moves a sixteenth of the way towards each answer it sees, in the encoder and the decoder alike, so they learn the chunk while it is being decoded. And there is not one probability per kind of question, but one per situation (the format calls these contexts): “same as west?” has a separate probability for every depth, combined with whether the west and south neighbours agree with each other, and whether the previous record agreed with either of them. Deep inside a city, where both neighbours and the previous record are all the same place, the “same as west?” probability for that situation is pushed to its limit, and a record costs about 0.005 bits. That is a run-length encoding without ever writing down a run length.

Why not Huffman? A Huffman code spends at least one whole bit for every symbol it writes. Huffman coding the depths and the palette indices separately costs at least two bits per record, already more than the 1.34 bits version 2 spends on both. Claude’s experiments showed that range coding gave the largest compression.

bound = (range >> 12) × p
if code < bound:  bit = 0,  range = bound
else:             bit = 1,  code = code − bound,  range = range − bound