From a50c43835d4cfb46ab39a4c70840c437fb99263a Mon Sep 17 00:00:00 2001 From: Brandon Liu Date: Sun, 24 Jul 2022 18:40:43 +0800 Subject: [PATCH] TileCoord supports up to zoom 15 using alternate ordering (#266) --- CONTRIBUTING.md | 6 + .../benchmarks/BenchmarkTileCoord.java | 38 ++++ .../onthegomap/planetiler/geo/TileCoord.java | 196 ++++++++++++++---- .../planetiler/mbtiles/MbtilesWriter.java | 5 +- .../planetiler/geo/TileCoordTest.java | 86 ++++---- 5 files changed, 244 insertions(+), 87 deletions(-) create mode 100644 planetiler-benchmarks/src/main/java/com/onthegomap/planetiler/benchmarks/BenchmarkTileCoord.java diff --git a/CONTRIBUTING.md b/CONTRIBUTING.md index fe9d6760..95f66159 100644 --- a/CONTRIBUTING.md +++ b/CONTRIBUTING.md @@ -23,6 +23,12 @@ To set up your local development environment: - on windows: `mvnw.cmd clean test` - or if you already have maven installed globally on your machine: `mvn clean test` - to run just one test e.g. `GeoUtilsTest`: `./mvnw -pl planetiler-core -Dtest=GeoUtilsTest test` + - to run benchmarks e.g. `BenchmarkTileCoord`: + +```sh +./scripts/build.sh +java -cp planetiler-dist/target/planetiler-dist-*-with-deps.jar com.onthegomap.planetiler.benchmarks.BenchmarkTileCoord +``` GitHub Workflows will run regression tests on any pull request. diff --git a/planetiler-benchmarks/src/main/java/com/onthegomap/planetiler/benchmarks/BenchmarkTileCoord.java b/planetiler-benchmarks/src/main/java/com/onthegomap/planetiler/benchmarks/BenchmarkTileCoord.java new file mode 100644 index 00000000..fffd005f --- /dev/null +++ b/planetiler-benchmarks/src/main/java/com/onthegomap/planetiler/benchmarks/BenchmarkTileCoord.java @@ -0,0 +1,38 @@ +package com.onthegomap.planetiler.benchmarks; + +import static io.prometheus.client.Collector.NANOSECONDS_PER_SECOND; + +import com.onthegomap.planetiler.geo.TileCoord; +import com.onthegomap.planetiler.stats.Timer; +import com.onthegomap.planetiler.util.Format; + +public class BenchmarkTileCoord { + + public static void main(String[] args) { + for (int i = 0; i < 3; i++) { + var timer = Timer.start(); + int num = 0; + for (int z = 0; z <= 14; z++) { + int max = 1 << z; + for (int x = 0; x < max; x++) { + for (int y = 0; y < max; y++) { + int encoded = TileCoord.encode(x, y, z); + int decoded = TileCoord.decode(encoded).encoded(); + // make sure we use the result so it doesn't get jit'ed-out + if (encoded != decoded) { + System.err.println("Error on " + z + "/" + x + "/" + y); + } + num++; + } + } + } + System.err.println( + "z0-z14 took " + + Format.defaultInstance().duration(timer.stop().elapsed().wall()) + " (" + + Format.defaultInstance() + .numeric(num * 1d / (timer.stop().elapsed().wall().toNanos() / NANOSECONDS_PER_SECOND)) + + "/s)" + ); + } + } +} diff --git a/planetiler-core/src/main/java/com/onthegomap/planetiler/geo/TileCoord.java b/planetiler-core/src/main/java/com/onthegomap/planetiler/geo/TileCoord.java index be209361..8c429244 100644 --- a/planetiler-core/src/main/java/com/onthegomap/planetiler/geo/TileCoord.java +++ b/planetiler-core/src/main/java/com/onthegomap/planetiler/geo/TileCoord.java @@ -1,6 +1,5 @@ package com.onthegomap.planetiler.geo; -import com.onthegomap.planetiler.mbtiles.Mbtiles; import com.onthegomap.planetiler.util.Format; import javax.annotation.concurrent.Immutable; import org.locationtech.jts.geom.Coordinate; @@ -9,29 +8,21 @@ import org.locationtech.jts.geom.CoordinateXY; /** * The coordinate of a slippy map tile. *

- * In order to encode into a 32-bit integer, only zoom levels {@code <= 14} are supported since we need 4 bits for the - * zoom-level, and 14 bits each for the x/y coordinates. + * Tile coords are sorted by consecutive Z levels in ascending order: 0 coords for z=0, 4 coords for z=1, etc. TMS + * order: tiles in a level are sorted by x ascending, y descending to match the ordering of the MBTiles sqlite index. + * Hilbert order: tiles in a level are ordered on the Hilbert curve with the first coordinate at the tip left. *

- * Tiles are ordered by z ascending, x ascending, y descending to match index ordering of {@link Mbtiles} sqlite - * database. * * @param encoded the tile ID encoded as a 32-bit integer * @param x x coordinate of the tile where 0 is the western-most tile just to the east the international date line * and 2^z-1 is the eastern-most tile * @param y y coordinate of the tile where 0 is the northern-most tile and 2^z-1 is the southern-most tile - * @param z zoom level ({@code <= 14}) + * @param z zoom level ({@code <= 15}) */ @Immutable public record TileCoord(int encoded, int x, int y, int z) implements Comparable { - // TODO: support higher than z14 - // z15 could theoretically fit into a 32-bit integer but needs a different packing strategy - // z16+ would need more space - // also need to remove hardcoded z14 limits - - private static final int XY_MASK = (1 << 14) - 1; - public TileCoord { - assert z <= 14; + assert z <= 15; } public static TileCoord ofXYZ(int x, int y, int z) { @@ -39,10 +30,19 @@ public record TileCoord(int encoded, int x, int y, int z) implements Comparable< } public static TileCoord decode(int encoded) { - int z = (encoded >> 28) + 8; - int x = (encoded >> 14) & XY_MASK; - int y = ((1 << z) - 1) - ((encoded) & XY_MASK); - return new TileCoord(encoded, x, y, z); + int acc = 0; + int tmpZ = 0; + while (true) { + int numTiles = (1 << tmpZ) * (1 << tmpZ); + if (acc + numTiles > encoded) { + int position = encoded - acc; + // long xy = hilbertPositionToXY(tmpZ, position); + long xy = tmsPositionToXY(tmpZ, position); + return new TileCoord(encoded, (int) (xy >>> 32 & 0xFFFFFFFFL), (int) (xy & 0xFFFFFFFFL), tmpZ); + } + acc += numTiles; + tmpZ++; + } } /** Returns the tile containing a latitude/longitude coordinate at a given zoom level. */ @@ -53,31 +53,13 @@ public record TileCoord(int encoded, int x, int y, int z) implements Comparable< return TileCoord.ofXYZ((int) Math.floor(x), (int) Math.floor(y), zoom); } - private static int encode(int x, int y, int z) { - int max = 1 << z; - if (x >= max) { - x %= max; + public static int encode(int x, int y, int z) { + int acc = 0; + for (int tmpZ = 0; tmpZ < z; tmpZ++) { + acc += (1 << tmpZ) * (1 << tmpZ); } - if (x < 0) { - x += max; - } - if (y < 0) { - y = 0; - } - if (y >= max) { - y = max - 1; - } - // since most significant bit is treated as the sign bit, make: - // z0-7 get encoded from 8 (0b1000) to 15 (0b1111) - // z8-14 get encoded from 0 (0b0000) to 6 (0b0110) - // so that encoded tile coordinates are ordered by zoom level - if (z < 8) { - z += 8; - } else { - z -= 8; - } - y = max - 1 - y; - return (z << 28) | (x << 14) | y; + // return acc + hilbertXYToPosition(z, x, y); + return acc + tmsXYToPosition(z, x, y); } @Override @@ -104,6 +86,19 @@ public record TileCoord(int encoded, int x, int y, int z) implements Comparable< return "{x=" + x + " y=" + y + " z=" + z + '}'; } + public double progressOnLevel() { + int acc = 0; + int tmpZ = 0; + while (true) { + int numTiles = (1 << tmpZ) * (1 << tmpZ); + if (acc + numTiles > encoded) { + return (encoded - acc) / (double) numTiles; + } + acc += numTiles; + tmpZ++; + } + } + @Override public int compareTo(TileCoord o) { return Long.compare(encoded, o.encoded); @@ -132,4 +127,121 @@ public record TileCoord(int encoded, int x, int y, int z) implements Comparable< double y = GeoUtils.getWorldY(lat) * factor; return new CoordinateXY((x - Math.floor(x)) * 256, (y - Math.floor(y)) * 256); } + + public static long tmsPositionToXY(int z, int pos) { + if (z == 0) + return 0; + int dim = 1 << z; + int x = pos / dim; + int y = dim - 1 - (pos % dim); + return ((long) x << 32) | y; + } + + public static int tmsXYToPosition(int z, int x, int y) { + int dim = 1 << z; + return x * dim + (dim - 1 - y); + } + + // hilbert implementation (not currently used) + // Fast Hilbert curve algorithm by http://threadlocalmutex.com/ + // Ported from C++ https://github.com/rawrunprotected/hilbert_curves (public domain) + private static int deinterleave(int tx) { + tx = tx & 0x55555555; + tx = (tx | (tx >>> 1)) & 0x33333333; + tx = (tx | (tx >>> 2)) & 0x0F0F0F0F; + tx = (tx | (tx >>> 4)) & 0x00FF00FF; + tx = (tx | (tx >>> 8)) & 0x0000FFFF; + return tx; + } + + private static int interleave(int tx) { + tx = (tx | (tx << 8)) & 0x00FF00FF; + tx = (tx | (tx << 4)) & 0x0F0F0F0F; + tx = (tx | (tx << 2)) & 0x33333333; + tx = (tx | (tx << 1)) & 0x55555555; + return tx; + } + + private static int prefixScan(int tx) { + tx = (tx >>> 8) ^ tx; + tx = (tx >>> 4) ^ tx; + tx = (tx >>> 2) ^ tx; + tx = (tx >>> 1) ^ tx; + return tx; + } + + private static long hilbertPositionToXY(int z, int pos) { + pos = pos << (32 - 2 * z); + + int i0 = deinterleave(pos); + int i1 = deinterleave(pos >>> 1); + + int t0 = (i0 | i1) ^ 0xFFFF; + int t1 = i0 & i1; + + int prefixT0 = prefixScan(t0); + int prefixT1 = prefixScan(t1); + + int a = (((i0 ^ 0xFFFF) & prefixT1) | (i0 & prefixT0)); + + int resultX = (a ^ i1) >>> (16 - z); + int resultY = (a ^ i0 ^ i1) >>> (16 - z); + return ((long) resultX << 32) | resultY; + } + + private static int hilbertXYToIndex(int z, int x, int y) { + x = x << (16 - z); + y = y << (16 - z); + + int hA, hB, hC, hD; + + int a1 = x ^ y; + int b1 = 0xFFFF ^ a1; + int c1 = 0xFFFF ^ (x | y); + int d1 = x & (y ^ 0xFFFF); + + hA = a1 | (b1 >>> 1); + hB = (a1 >>> 1) ^ a1; + + hC = ((c1 >>> 1) ^ (b1 & (d1 >>> 1))) ^ c1; + hD = ((a1 & (c1 >>> 1)) ^ (d1 >>> 1)) ^ d1; + + int a2 = hA; + int b2 = hB; + int c2 = hC; + int d2 = hD; + + hA = ((a2 & (a2 >>> 2)) ^ (b2 & (b2 >>> 2))); + hB = ((a2 & (b2 >>> 2)) ^ (b2 & ((a2 ^ b2) >>> 2))); + + hC ^= ((a2 & (c2 >>> 2)) ^ (b2 & (d2 >>> 2))); + hD ^= ((b2 & (c2 >>> 2)) ^ ((a2 ^ b2) & (d2 >>> 2))); + + int a3 = hA; + int b3 = hB; + int c3 = hC; + int d3 = hD; + + hA = ((a3 & (a3 >>> 4)) ^ (b3 & (b3 >>> 4))); + hB = ((a3 & (b3 >>> 4)) ^ (b3 & ((a3 ^ b3) >>> 4))); + + hC ^= ((a3 & (c3 >>> 4)) ^ (b3 & (d3 >>> 4))); + hD ^= ((b3 & (c3 >>> 4)) ^ ((a3 ^ b3) & (d3 >>> 4))); + + int a4 = hA; + int b4 = hB; + int c4 = hC; + int d4 = hD; + + hC ^= ((a4 & (c4 >>> 8)) ^ (b4 & (d4 >>> 8))); + hD ^= ((b4 & (c4 >>> 8)) ^ ((a4 ^ b4) & (d4 >>> 8))); + + int a = hC ^ (hC >>> 1); + int b = hD ^ (hD >>> 1); + + int i0 = x ^ y; + int i1 = b | (0xFFFF ^ (i0 | a)); + + return ((interleave(i1) << 1) | interleave(i0)) >>> (32 - 2 * z); + } } diff --git a/planetiler-core/src/main/java/com/onthegomap/planetiler/mbtiles/MbtilesWriter.java b/planetiler-core/src/main/java/com/onthegomap/planetiler/mbtiles/MbtilesWriter.java index 31b52e40..3b14db18 100644 --- a/planetiler-core/src/main/java/com/onthegomap/planetiler/mbtiles/MbtilesWriter.java +++ b/planetiler-core/src/main/java/com/onthegomap/planetiler/mbtiles/MbtilesWriter.java @@ -203,12 +203,9 @@ public class MbtilesWriter { if (lastTile == null) { blurb = "n/a"; } else { - var extentForZoom = config.bounds().tileExtents().getForZoom(lastTile.z()); - int zMinX = extentForZoom.minX(); - int zMaxX = extentForZoom.maxX(); blurb = "%d/%d/%d (z%d %s%%) %s".formatted( lastTile.z(), lastTile.x(), lastTile.y(), - lastTile.z(), (100 * (lastTile.x() + 1 - zMinX)) / (zMaxX - zMinX), + lastTile.z(), 100 * lastTile.progressOnLevel(), lastTile.getDebugUrl() ); } diff --git a/planetiler-core/src/test/java/com/onthegomap/planetiler/geo/TileCoordTest.java b/planetiler-core/src/test/java/com/onthegomap/planetiler/geo/TileCoordTest.java index ce7e4964..0adce294 100644 --- a/planetiler-core/src/test/java/com/onthegomap/planetiler/geo/TileCoordTest.java +++ b/planetiler-core/src/test/java/com/onthegomap/planetiler/geo/TileCoordTest.java @@ -1,57 +1,61 @@ package com.onthegomap.planetiler.geo; import static org.junit.jupiter.api.Assertions.assertEquals; -import static org.junit.jupiter.api.Assertions.assertThrows; -import static org.junit.jupiter.api.Assertions.fail; -import org.junit.jupiter.api.Test; import org.junit.jupiter.params.ParameterizedTest; import org.junit.jupiter.params.provider.CsvSource; + class TileCoordTest { @ParameterizedTest @CsvSource({ - "0,0,0", - "0,0,1", - "0,1,1", - "1,1,1", - "100,100,14", + "0,0,0,0", + "0,1,1,1", + "0,0,1,2", + "1,1,1,3", + "1,0,1,4", + "0,3,2,5", + "0,2,2,6", + "0,1,2,7", + "0,0,2,8", + "1,3,2,9", + "1,2,2,10", + "1,1,2,11", + "1,0,2,12", + "2,3,2,13", + "2,2,2,14", + "2,1,2,15", + "2,0,2,16", + "3,3,2,17", + "3,2,2,18", + "3,1,2,19", + "3,0,2,20", + "0,0,15,357946708", + "0,32767,15,357913941", + "32767,0,15,1431655764", + "32767,32767,15,1431622997" }) - void testTileCoord(int x, int y, int z) { - TileCoord coord1 = TileCoord.ofXYZ(x, y, z); - TileCoord coord2 = TileCoord.decode(coord1.encoded()); - assertEquals(coord1.x(), coord2.x(), "x"); - assertEquals(coord1.y(), coord2.y(), "y"); - assertEquals(coord1.z(), coord2.z(), "z"); - assertEquals(coord1, coord2); + void testTileOrder(int x, int y, int z, int i) { + int encoded = TileCoord.ofXYZ(x, y, z).encoded(); + assertEquals(i, encoded); + TileCoord decoded = TileCoord.decode(i); + assertEquals(decoded.x(), x, "x"); + assertEquals(decoded.y(), y, "y"); + assertEquals(decoded.z(), z, "z"); } - @Test - void testTileSortOrderRespectZ() { - int last = Integer.MIN_VALUE; - for (int z = 0; z <= 14; z++) { - int encoded = TileCoord.ofXYZ(0, 0, z).encoded(); - if (encoded < last) { - fail("encoded value for z" + (z - 1) + " (" + last + ") is not less than z" + z + " (" + encoded + ")"); - } - last = encoded; - } - } - - @Test - void testTileSortOrderFlipY() { - for (int z = 1; z <= 14; z++) { - int encoded1 = TileCoord.ofXYZ(0, 1, z).encoded(); - int encoded2 = TileCoord.ofXYZ(0, 0, z).encoded(); - if (encoded2 < encoded1) { - fail("encoded value for y=1 is not less than y=0 at z=" + z); - } - } - } - - @Test - void testThrowsPastZ14() { - assertThrows(AssertionError.class, () -> TileCoord.ofXYZ(0, 0, 15)); + @ParameterizedTest + @CsvSource({ + "0,0,0,0", + "0,1,1,0", + "0,0,1,0.25", + "1,1,1,0.5", + "1,0,1,0.75", + "0,3,2,0" + }) + void testTileProgressOnLevel(int x, int y, int z, double p) { + double progress = TileCoord.ofXYZ(x, y, z).progressOnLevel(); + assertEquals(p, progress); } }