diff --git a/clip.cpp b/clip.cpp index d4a6f9ce..9b778391 100644 --- a/clip.cpp +++ b/clip.cpp @@ -571,8 +571,7 @@ drawvec remove_noop(drawvec geom, int type, int shift) { } if (geom[i + 1].op == VT_CLOSEPATH) { - // followed by closepath: not possible - fprintf(stderr, "Shouldn't happen\n"); + // followed by closepath: only possible after close_poly() i++; // also remove unused closepath continue; } @@ -2469,3 +2468,87 @@ drawvec fix_polygon(const drawvec &geom, bool use_winding, bool reverse_winding) return out; } + +bool line_is_too_small(drawvec const &geometry, int z, int detail) { + if (geometry.size() == 0) { + return true; + } + + long long x = 0, y = 0; + for (auto &g : geometry) { + if (g.op == VT_MOVETO) { + x = std::llround((double) g.x / (1LL << (32 - detail - z))); + y = std::llround((double) g.y / (1LL << (32 - detail - z))); + } else { + long long xx = std::llround((double) g.x / (1LL << (32 - detail - z))); + long long yy = std::llround((double) g.y / (1LL << (32 - detail - z))); + + if (xx != x || yy != y) { + return false; + } + } + } + + return true; +} + +void coalesce_polygon(drawvec &geom, bool scale_up) { + // wagyu should be able to straightforwardly handle + // anything under a few hundred thousand vertices + if (geom.size() < 100000) { + geom = clean_or_clip_poly(geom, 0, 0, false, scale_up); + return; + } + + // These geometries were assembled in geometric order, + // so sub-batches of them should hopefully union into + // reasonable sets. + // + // Find the first outer ring after halfway point. + + for (size_t i = geom.size() / 2; i < geom.size(); i++) { + if (geom[i].op == VT_MOVETO) { + size_t j; + for (j = i + 1; j < geom.size(); j++) { + if (geom[j].op != VT_LINETO) { + break; + } + } + + if (get_area(geom, i, j) > 0) { + // If we have an outer ring, split there + // and coalesce the two halves + + // Copy second half to new vector + std::vector geom2; + geom2.resize(geom.size() - i); + for (size_t k = i; k < geom.size(); k++) { + geom2[k - i] = geom[k]; + } + + // Resize vector to include only first half + geom.resize(i); + + // Clean each half individually + coalesce_polygon(geom, scale_up); + coalesce_polygon(geom2, scale_up); + + // Copy second half back with first + size_t brk = geom.size(); + geom.resize(brk + geom2.size()); + for (size_t k = 0; k < geom2.size(); k++) { + geom[brk + k] = geom2[k]; + } + geom2.clear(); + + // Clean the combined geometry + geom = clean_or_clip_poly(geom, 0, 0, false, scale_up); + } + + i = j - 1; + } + } + + // Can't find a breakpoint; take what we can get. + geom = clean_or_clip_poly(geom, 0, 0, false, scale_up); +} diff --git a/geometry.hpp b/geometry.hpp index 73a621ba..454cea61 100644 --- a/geometry.hpp +++ b/geometry.hpp @@ -173,4 +173,7 @@ void get_quadkey_bounds(long long xmin, long long ymin, long long xmax, long lon clipbbox parse_clip_poly(std::string arg); +bool line_is_too_small(drawvec const &geometry, int z, int detail); +void coalesce_polygon(drawvec &geom, bool scale_up); + #endif diff --git a/tile.cpp b/tile.cpp index dc4b8d7d..b7af6ba2 100644 --- a/tile.cpp +++ b/tile.cpp @@ -639,7 +639,7 @@ static double simplify_feature(serial_feature *p, drawvec const &shared_nodes, n // unioned exactly // // don't try to scale up because these are still world coordinates - geom = clean_or_clip_poly(geom, 0, 0, false, false); + coalesce_polygon(geom, false); } // continues to simplify to line_detail even if we have extra detail @@ -695,7 +695,7 @@ static void *simplification_worker(void *v) { if (!a->trying_to_stop_early) { // we can try scaling up because this is now tile scale - geom = clean_or_clip_poly(geom, 0, 0, false, true); + coalesce_polygon(geom, true); if (additional[A_DEBUG_POLYGON]) { check_polygon(geom); } @@ -1551,26 +1551,6 @@ bool find_feature_to_accumulate_onto(std::vector return false; } -static bool line_is_too_small(drawvec const &geometry, int z, int detail) { - if (geometry.size() == 0) { - return true; - } - - long long x = std::round((double) geometry[0].x / (1LL << (32 - detail - z))); - long long y = std::round((double) geometry[0].y / (1LL << (32 - detail - z))); - - for (auto &g : geometry) { - long long xx = std::round((double) g.x / (1LL << (32 - detail - z))); - long long yy = std::round((double) g.y / (1LL << (32 - detail - z))); - - if (xx != x || yy != y) { - return false; - } - } - - return true; -} - // Keep only a sample of 100K extents for feature dropping, // to avoid spending lots of memory on a complete list when there are // hundreds of millions of features. @@ -2201,7 +2181,7 @@ long long write_tile(decompressor *geoms, std::atomic *geompos_in, ch drawvec to_clean = features[simplified_geometry_through]->geometry; // don't scale up because this is still world coordinates - to_clean = clean_or_clip_poly(to_clean, 0, 0, false, false); + coalesce_polygon(to_clean, false); features[simplified_geometry_through]->geometry = std::move(to_clean); } } @@ -2470,7 +2450,7 @@ long long write_tile(decompressor *geoms, std::atomic *geompos_in, ch if (layer_features[x]->t == VT_POLYGON) { if (layer_features[x]->coalesced) { // we can try scaling up because this is tile coordinates - layer_features[x]->geometry = clean_or_clip_poly(layer_features[x]->geometry, 0, 0, false, true); + coalesce_polygon(layer_features[x]->geometry, true); } layer_features[x]->geometry = close_poly(layer_features[x]->geometry); diff --git a/unit.cpp b/unit.cpp index c8f36491..5e0bfa5d 100644 --- a/unit.cpp +++ b/unit.cpp @@ -152,3 +152,12 @@ TEST_CASE("mvt_geometry bbox") { REQUIRE(start == 0x1c84fc0000000000); REQUIRE(end == 0x1c84ffffffffffff); } + +TEST_CASE("line_is_too_small") { + drawvec dv; + dv.emplace_back(VT_MOVETO, 4243099709, 2683872952); + dv.emplace_back(VT_LINETO, 4243102487, 2683873977); + dv.emplace_back(VT_MOVETO, -51867587, 2683872952); + dv.emplace_back(VT_LINETO, -51864809, 2683873977); + REQUIRE(line_is_too_small(dv, 0, 10)); +}