1#pragma once
2
3#include <mapbox/geojsonvt/convert.hpp>
4#include <mapbox/geojsonvt/tile.hpp>
5#include <mapbox/geojsonvt/types.hpp>
6#include <mapbox/geojsonvt/wrap.hpp>
7
8#include <chrono>
9#include <cmath>
10#include <map>
11#include <unordered_map>
12
13namespace mapbox {
14namespace geojsonvt {
15
16using geometry = mapbox::geometry::geometry<double>;
17using feature = mapbox::geometry::feature<double>;
18using feature_collection = mapbox::geometry::feature_collection<double>;
19using geometry_collection = mapbox::geometry::geometry_collection<double>;
20using geojson = mapbox::util::variant<geometry, feature, feature_collection>;
21
22struct ToFeatureCollection {
23 feature_collection operator()(const feature_collection& value) const {
24 return value;
25 }
26 feature_collection operator()(const feature& value) const {
27 return { value };
28 }
29 feature_collection operator()(const geometry& value) const {
30 return { { value } };
31 }
32};
33
34struct TileOptions {
35 // simplification tolerance (higher means simpler)
36 double tolerance = 3;
37
38 // tile extent
39 uint16_t extent = 4096;
40
41 // tile buffer on each side
42 uint16_t buffer = 64;
43};
44
45struct Options : TileOptions {
46 // max zoom to preserve detail on
47 uint8_t maxZoom = 18;
48
49 // max zoom in the tile index
50 uint8_t indexMaxZoom = 5;
51
52 // max number of points per tile in the tile index
53 uint32_t indexMaxPoints = 100000;
54};
55
56const Tile empty_tile{};
57
58inline uint64_t toID(uint8_t z, uint32_t x, uint32_t y) {
59 return (((1ull << z) * y + x) * 32) + z;
60}
61
62inline const Tile geoJSONToTile(const geojson& geojson_,
63 uint8_t z,
64 uint32_t x,
65 uint32_t y,
66 const TileOptions& options = TileOptions(),
67 bool wrap = false,
68 bool clip = false) {
69
70 const auto features_ = geojson::visit(v: geojson_, f: ToFeatureCollection{});
71 auto z2 = 1u << z;
72 auto tolerance = (options.tolerance / options.extent) / z2;
73 auto features = detail::convert(features: features_, tolerance);
74 if (wrap) {
75 features = detail::wrap(features, buffer: double(options.buffer) / options.extent);
76 }
77 if (clip) {
78 const double p = options.buffer / options.extent;
79
80 const auto left = detail::clip<0>(features, k1: (x - p) / z2, k2: (x + 1 + p) / z2, minAll: -1, maxAll: 2);
81 features = detail::clip<1>(features: left, k1: (y - p) / z2, k2: (y + 1 + p) / z2, minAll: -1, maxAll: 2);
82 }
83 return detail::InternalTile({ features, z, x, y, options.extent, tolerance }).tile;
84}
85
86class GeoJSONVT {
87public:
88 const Options options;
89
90 GeoJSONVT(const mapbox::geometry::feature_collection<double>& features_,
91 const Options& options_ = Options())
92 : options(options_) {
93
94 const uint32_t z2 = 1u << options.maxZoom;
95
96 auto converted = detail::convert(features: features_, tolerance: (options.tolerance / options.extent) / z2);
97 auto features = detail::wrap(features: converted, buffer: double(options.buffer) / options.extent);
98
99 splitTile(features, z: 0, x: 0, y: 0);
100 }
101
102 GeoJSONVT(const geojson& geojson_, const Options& options_ = Options())
103 : GeoJSONVT(geojson::visit(v: geojson_, f: ToFeatureCollection{}), options_) {
104 }
105
106 std::map<uint8_t, uint32_t> stats;
107 uint32_t total = 0;
108
109 const Tile& getTile(const uint8_t z, const uint32_t x_, const uint32_t y) {
110
111 if (z > options.maxZoom)
112 throw std::runtime_error("Requested zoom higher than maxZoom: " + std::to_string(val: z));
113
114 const uint32_t z2 = 1u << z;
115 const uint32_t x = ((x_ % z2) + z2) % z2; // wrap tile x coordinate
116 const uint64_t id = toID(z, x, y);
117
118 auto it = tiles.find(x: id);
119 if (it != tiles.end())
120 return it->second.tile;
121
122 it = findParent(z, x, y);
123
124 if (it == tiles.end())
125 throw std::runtime_error("Parent tile not found");
126
127 // if we found a parent tile containing the original geometry, we can drill down from it
128 const auto& parent = it->second;
129
130 // drill down parent tile up to the requested one
131 splitTile(features: parent.source_features, z: parent.z, x: parent.x, y: parent.y, cz: z, cx: x, cy: y);
132
133 it = tiles.find(x: id);
134 if (it != tiles.end())
135 return it->second.tile;
136
137 it = findParent(z, x, y);
138 if (it == tiles.end())
139 throw std::runtime_error("Parent tile not found");
140
141 return empty_tile;
142 }
143
144 const std::unordered_map<uint64_t, detail::InternalTile>& getInternalTiles() const {
145 return tiles;
146 }
147
148private:
149 std::unordered_map<uint64_t, detail::InternalTile> tiles;
150
151 std::unordered_map<uint64_t, detail::InternalTile>::iterator
152 findParent(const uint8_t z, const uint32_t x, const uint32_t y) {
153 uint8_t z0 = z;
154 uint32_t x0 = x;
155 uint32_t y0 = y;
156
157 const auto end = tiles.end();
158 auto parent = end;
159
160 while ((parent == end) && (z0 != 0)) {
161 z0--;
162 x0 = x0 / 2;
163 y0 = y0 / 2;
164 parent = tiles.find(x: toID(z: z0, x: x0, y: y0));
165 }
166
167 return parent;
168 }
169
170 void splitTile(const detail::vt_features& features,
171 const uint8_t z,
172 const uint32_t x,
173 const uint32_t y,
174 const uint8_t cz = 0,
175 const uint32_t cx = 0,
176 const uint32_t cy = 0) {
177
178 const double z2 = 1u << z;
179 const uint64_t id = toID(z, x, y);
180
181 auto it = tiles.find(x: id);
182
183 if (it == tiles.end()) {
184 const double tolerance =
185 (z == options.maxZoom ? 0 : options.tolerance / (z2 * options.extent));
186
187 it = tiles
188 .emplace(args: id,
189 args: detail::InternalTile{ features, z, x, y, options.extent, tolerance })
190 .first;
191 stats[z] = (stats.count(x: z) ? stats[z] + 1 : 1);
192 total++;
193 // printf("tile z%i-%i-%i\n", z, x, y);
194 }
195
196 auto& tile = it->second;
197
198 if (features.empty())
199 return;
200
201 // if it's the first-pass tiling
202 if (cz == 0u) {
203 // stop tiling if we reached max zoom, or if the tile is too simple
204 if (z == options.indexMaxZoom || tile.tile.num_points <= options.indexMaxPoints) {
205 tile.source_features = features;
206 return;
207 }
208
209 } else { // drilldown to a specific tile;
210 // stop tiling if we reached base zoom
211 if (z == options.maxZoom)
212 return;
213
214 // stop tiling if it's our target tile zoom
215 if (z == cz) {
216 tile.source_features = features;
217 return;
218 }
219
220 // stop tiling if it's not an ancestor of the target tile
221 const double m = 1u << (cz - z);
222 if (x != static_cast<uint32_t>(std::floor(x: cx / m)) ||
223 y != static_cast<uint32_t>(std::floor(x: cy / m))) {
224 tile.source_features = features;
225 return;
226 }
227 }
228
229 const double p = 0.5 * options.buffer / options.extent;
230 const auto& min = tile.bbox.min;
231 const auto& max = tile.bbox.max;
232
233 const auto left = detail::clip<0>(features, k1: (x - p) / z2, k2: (x + 0.5 + p) / z2, minAll: min.x, maxAll: max.x);
234
235 splitTile(features: detail::clip<1>(features: left, k1: (y - p) / z2, k2: (y + 0.5 + p) / z2, minAll: min.y, maxAll: max.y), z: z + 1,
236 x: x * 2, y: y * 2, cz, cx, cy);
237 splitTile(features: detail::clip<1>(features: left, k1: (y + 0.5 - p) / z2, k2: (y + 1 + p) / z2, minAll: min.y, maxAll: max.y), z: z + 1,
238 x: x * 2, y: y * 2 + 1, cz, cx, cy);
239
240 const auto right =
241 detail::clip<0>(features, k1: (x + 0.5 - p) / z2, k2: (x + 1 + p) / z2, minAll: min.x, maxAll: max.x);
242
243 splitTile(features: detail::clip<1>(features: right, k1: (y - p) / z2, k2: (y + 0.5 + p) / z2, minAll: min.y, maxAll: max.y), z: z + 1,
244 x: x * 2 + 1, y: y * 2, cz, cx, cy);
245 splitTile(features: detail::clip<1>(features: right, k1: (y + 0.5 - p) / z2, k2: (y + 1 + p) / z2, minAll: min.y, maxAll: max.y), z: z + 1,
246 x: x * 2 + 1, y: y * 2 + 1, cz, cx, cy);
247
248 // if we sliced further down, no need to keep source geometry
249 tile.source_features = {};
250 }
251};
252
253} // namespace geojsonvt
254} // namespace mapbox
255

source code of qtlocation/src/3rdparty/mapbox-gl-native/deps/geojsonvt/6.5.1/include/mapbox/geojsonvt.hpp