Skip to content

Commit

Permalink
Simplify and optimize Hilbert tile ID <-> XYZ conversion (#527)
Browse files Browse the repository at this point in the history 
* Optimize Hilbert tile ID <-> XYZ conversion

* biome

* tile_id: more improvements

* tile_id(cpp): uint32
  • Loading branch information
@ciscorn
ciscorn authored Feb 21, 2025
1 parent caa8f1f commit 9b69ab4
Show file tree
Hide file tree
Showing 4 changed files with 122 additions and 157 deletions.
74 changes: 35 additions & 39 deletions cpp/pmtiles.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -321,35 +321,18 @@ uint64_t decode_varint(const char **data, const char *end) {
return decode_varint_impl(data, end);
}

void rotate(int64_t n, int64_t &x, int64_t &y, int64_t rx, int64_t ry) {
void rotate(int64_t n, uint32_t &x, uint32_t &y, uint32_t rx, uint32_t ry) {
if (ry == 0) {
if (rx == 1) {
if (rx != 0) {
x = n - 1 - x;
y = n - 1 - y;
}
int64_t t = x;
uint32_t t = x;
x = y;
y = t;
}
}

zxy t_on_level(uint8_t z, uint64_t pos) {
int64_t n = 1LL << z;
int64_t rx, ry, s, t = pos;
int64_t tx = 0;
int64_t ty = 0;

for (s = 1; s < n; s *= 2) {
rx = 1LL & (t / 2);
ry = 1LL & (t ^ rx);
rotate(s, tx, ty, rx, ry);
tx += s * rx;
ty += s * ry;
t /= 4;
}
return zxy(z, static_cast<int>(tx), static_cast<int>(ty));
}

int write_varint(std::back_insert_iterator<std::string> data, uint64_t value) {
int n = 1;

Expand Down Expand Up @@ -405,16 +388,31 @@ entryv3 find_tile(const std::vector<entryv3> &entries, uint64_t tile_id) {

} // end anonymous namespace

// Note: std::bit_width is available in C++20 and later.
static uint8_t bit_width(uint64_t n) {
uint8_t count = 0;
while (n > 0) { count++; n >>= 1; }
return count;
}

inline zxy tileid_to_zxy(uint64_t tileid) {
uint64_t acc = 0;
for (uint8_t t_z = 0; t_z < 32; t_z++) {
uint64_t num_tiles = (1LL << t_z) * (1LL << t_z);
if (acc + num_tiles > tileid) {
return t_on_level(t_z, tileid - acc);
}
acc += num_tiles;
if (tileid > 6148914691236517204) {
throw std::overflow_error("tile zoom exceeds 64-bit limit");
}
throw std::overflow_error("tile zoom exceeds 64-bit limit");
uint8_t z = (bit_width(3 * tileid + 1) - 1) / 2;
uint64_t acc = ((1L << (z * 2)) - 1) / 3;
uint64_t pos = tileid - acc;
uint32_t x = 0, y = 0;
for (uint8_t a = 0; a < z; a++) {
uint64_t s = 1 << a;
uint32_t rx = s & (pos / 2);
uint32_t ry = s & (pos ^ rx);
rotate(s, x, y, rx, ry);
pos >>= 1;
x += rx;
y += ry;
}
return zxy(z, static_cast<int>(x), static_cast<int>(y));
}

inline uint64_t zxy_to_tileid(uint8_t z, uint32_t x, uint32_t y) {
Expand All @@ -424,19 +422,17 @@ inline uint64_t zxy_to_tileid(uint8_t z, uint32_t x, uint32_t y) {
if (x > (1U << z) - 1U || y > (1U << z) - 1U) {
throw std::overflow_error("tile x/y outside zoom level bounds");
}
uint64_t acc = 0;
for (uint8_t t_z = 0; t_z < z; t_z++) acc += (1LL << t_z) * (1LL << t_z);
int64_t n = 1LL << z;
int64_t rx, ry, s, d = 0;
int64_t tx = x;
int64_t ty = y;
for (s = n / 2; s > 0; s /= 2) {
rx = (tx & s) > 0;
ry = (ty & s) > 0;
d += s * s * ((3LL * rx) ^ ry);
uint64_t acc = ((1LL << (z * 2U)) - 1) / 3;
uint32_t tx = x, ty = y;
int a = z - 1;
for (uint32_t s = 1LL << a; s > 0; s >>= 1) {
uint32_t rx = s & tx;
uint32_t ry = s & ty;
rotate(s, tx, ty, rx, ry);
acc += ((3LL * rx) ^ ry) << a;
a--;
}
return acc + d;
return acc;
}

// returns an uncompressed byte buffer
Expand Down
3 changes: 0 additions & 3 deletions cpp/test.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -72,7 +72,6 @@ MU_TEST(test_roundtrip) {
} catch (const std::runtime_error &e) {
caught = true;
}

mu_check(caught);

caught = false;
Expand All @@ -81,7 +80,6 @@ MU_TEST(test_roundtrip) {
} catch (const std::runtime_error &e) {
caught = true;
}

mu_check(caught);

caught = false;
Expand All @@ -90,7 +88,6 @@ MU_TEST(test_roundtrip) {
} catch (const std::runtime_error &e) {
caught = true;
}

mu_check(caught);

caught = false;
Expand Down
105 changes: 46 additions & 59 deletions js/src/index.ts
View file
Original file line number Diff line number Diff line change
Expand Up @@ -55,88 +55,75 @@ export function readVarint(p: BufferPosition): number {
return readVarintRemainder(val, p);
}

function rotate(n: number, xy: number[], rx: number, ry: number): void {
function rotate(
n: number,
x: number,
y: number,
rx: number,
ry: number
): [number, number] {
if (ry === 0) {
if (rx === 1) {
xy[0] = n - 1 - xy[0];
xy[1] = n - 1 - xy[1];
if (rx !== 0) {
return [n - 1 - y, n - 1 - x];
}
const t = xy[0];
xy[0] = xy[1];
xy[1] = t;
return [y, x];
}
return [x, y];
}

function idOnLevel(z: number, pos: number): [number, number, number] {
const n = 2 ** z;
let rx = pos;
let ry = pos;
let t = pos;
const xy = [0, 0];
let s = 1;
while (s < n) {
rx = 1 & (t / 2);
ry = 1 & (t ^ rx);
rotate(s, xy, rx, ry);
xy[0] += s * rx;
xy[1] += s * ry;
t = t / 4;
s *= 2;
}
return [z, xy[0], xy[1]];
}

const tzValues: number[] = [
0, 1, 5, 21, 85, 341, 1365, 5461, 21845, 87381, 349525, 1398101, 5592405,
22369621, 89478485, 357913941, 1431655765, 5726623061, 22906492245,
91625968981, 366503875925, 1466015503701, 5864062014805, 23456248059221,
93824992236885, 375299968947541, 1501199875790165,
];

/**
* Convert Z,X,Y to a Hilbert TileID.
*/
export function zxyToTileId(z: number, x: number, y: number): number {
if (z > 26) {
throw new Error("Tile zoom level exceeds max safe number limit (26)");
}
if (x > 2 ** z - 1 || y > 2 ** z - 1) {
if (x >= 1 << z || y >= 1 << z) {
throw new Error("tile x/y outside zoom level bounds");
}
let acc = ((1 << z) * (1 << z) - 1) / 3;
let a = z - 1;
let [tx, ty] = [x, y];
for (let s = 1 << a; s > 0; s >>= 1) {
const rx = tx & s;
const ry = ty & s;
acc += ((3 * rx) ^ ry) * (1 << a);
[tx, ty] = rotate(s, tx, ty, rx, ry);
a--;
}
return acc;
}

const acc = tzValues[z];
const n = 2 ** z;
let rx = 0;
let ry = 0;
let d = 0;
const xy = [x, y];
let s = n / 2;
while (s > 0) {
rx = (xy[0] & s) > 0 ? 1 : 0;
ry = (xy[1] & s) > 0 ? 1 : 0;
d += s * s * ((3 * rx) ^ ry);
rotate(s, xy, rx, ry);
s = s / 2;
function tileIdToZ(i: number): number {
const c = 3 * i + 1;
if (c < 0x100000000) {
return 31 - Math.clz32(c);
}
return acc + d;
return 63 - Math.clz32(c / 0x100000000);
}

/**
* Convert a Hilbert TileID to Z,X,Y.
*/
export function tileIdToZxy(i: number): [number, number, number] {
let acc = 0;
const z = 0;

for (let z = 0; z < 27; z++) {
const numTiles = (0x1 << z) * (0x1 << z);
if (acc + numTiles > i) {
return idOnLevel(z, i - acc);
}
acc += numTiles;
const z = tileIdToZ(i) >> 1;
if (z > 26)
throw new Error("Tile zoom level exceeds max safe number limit (26)");
const acc = ((1 << z) * (1 << z) - 1) / 3;

let t = i - acc;
let x = 0;
let y = 0;
const n = 1 << z;
for (let s = 1; s < n; s <<= 1) {
const rx = s & (t / 2);
const ry = s & (t ^ rx);
[x, y] = rotate(s, x, y, rx, ry);
t = t / 2;
x += rx;
y += ry;
}

throw new Error("Tile zoom level exceeds max safe number limit (26)");
return [z, x, y];
}

/**
Expand Down
Loading

0 comments on commit 9b69ab4

Please # to comment.