diff --git a/CMakeLists.txt b/CMakeLists.txt index 35608283..46ffd341 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -136,6 +136,7 @@ file(GLOB tilemaker_src_files src/config_validator.cpp src/coordinates.cpp src/coordinates_geom.cpp + src/declutter.cpp src/external/streamvbyte_decode.c src/external/streamvbyte_encode.c src/external/streamvbyte_zigzag.c diff --git a/Makefile b/Makefile index 298881bf..6cbf6580 100644 --- a/Makefile +++ b/Makefile @@ -100,6 +100,7 @@ tilemaker: \ src/config_validator.o \ src/coordinates_geom.o \ src/coordinates.o \ + src/declutter.o \ src/external/streamvbyte_decode.o \ src/external/streamvbyte_encode.o \ src/external/streamvbyte_zigzag.o \ diff --git a/docs/CONFIGURATION.md b/docs/CONFIGURATION.md index 348eacea..82fad41b 100644 --- a/docs/CONFIGURATION.md +++ b/docs/CONFIGURATION.md @@ -86,6 +86,9 @@ You can add optional parameters to layers: * `combine_lines_below` - whether to merge all linestrings in the tile with the same attributes. If not defined, the global setting `combine_below` will be used. * `combine_points` - merge points with the same attributes (defaults to `true`: specify `false` to disable) * `z_order_ascending` - sort features in ascending order by a numeric value set in the Lua processing script (defaults to `true`: specify `false` for descending order) +* `declutter_below` - thin out point features below this zoom level, so only the most important appear when zoomed out (see 'Decluttering point features' below) +* `declutter_distance` - how far apart (in pixels) decluttered features should be kept - defaults to 40 +* `declutter_threshold` - the score a feature needs to appear at the layer's `minzoom`; this halves at each subsequent zoom level - defaults to 0 (no threshold) `write_to` enables you to combine different layer specs within one outputted layer. For example: @@ -98,6 +101,30 @@ This would combine the `roads` (z12-14) and `low_roads` (z9-11) layers into a si (See also 'Shapefiles and GeoJSON' below.) +### Decluttering point features + +For more effective maps at lower zoom levels, tilemaker allows you to prioritise the most important point features - such as the largest cities or the highest peaks. + +In your Lua profile, give each point a `Score()`, typically calculated from OSM tags such as `population` or `ele`. Then set three `declutter` properties on the layer (see above) to determine which zoom levels this applies at, and the distance between each selected feature. + + "place": { + "minzoom": 4, "maxzoom": 14, + "declutter_below": 11, "declutter_distance": 40, "declutter_threshold": 50000 + } + +```lua + local place = Find("place") + if place=="city" or place=="town" or place=="village" then + Layer("place", false) + Attribute("name", Find("name")) + Score(tonumber(Find("population")) or 0) + end +``` + +Starting at the layer's `minzoom`, tilemaker takes the features whose score is at least `declutter_threshold` (highest score first), and places each one that isn't within `declutter_distance` of a feature already placed. It then moves up a zoom level, halving the threshold so that less important features are included. Anything still unplaced by the time we reach `declutter_below` is written from that zoom level upwards. + +A feature whose `MinZoom()` is above the zoom being considered isn't eligible for it. Only point features are decluttered, not lines or polygons. + ### Additional metadata Tilemaker writes a `json` metadata field containing a `vector_layers` key, whose value is an array of JSON objects describing each layer and its attributes. This is part of the MBTiles 1.3 spec and required by certain clients. @@ -171,6 +198,7 @@ To do that, you use these methods: * `IsMultiPolygon()`: returns true if the current object is a multipolygon. * `ZOrder(number)`: Set a numeric value (default 0) used to sort features within a layer. Use this feature to ensure a proper rendering order if the rendering engine itself does not support sorting. Sorting is not supported across layers merged with `write_to`. Features with different z-order are not merged if `combine_below`, `combine_lines_below` or `combine_polygons_below` is used. Use this in conjunction with `feature_limit` to only write the most important (highest z-order) features within a tile. (Values can be -50,000,000 to 50,000,000 and are lossy, particularly beyond -1000 to 1000.) * `MinZoom(zoom)`: set the minimum zoom level (0-15) at which this object will be written. Note that the JSON layer configuration minimum still applies (so `:MinZoom(5)` will have no effect if your layer only starts at z6). +* `Score(number)`: set how important this point feature is (a whole number; anything else is rounded down), so that tilemaker can decide which zoom level to show it from. Only meaningful for points (i.e. features written with `Layer` for a node, or `LayerAsCentroid`) in a layer with `declutter_below` set - see 'Decluttering point features' below. * `Length()` and `Area()`: return the length (metres)/area (square metres) of the current object. Requires Boost 1.67+. * `Centroid()`: return the lat/lon of the centre of the current object as a two-element Lua table (element 1 is lat, 2 is lon). @@ -233,6 +261,8 @@ Limited Lua transformations are available for these files. You can supply an `at To set the minimum zoom level at which an individual feature is rendered, use `attribute_function` to set a `_minzoom` value in your return table. +Similarly, to declutter shapefile/GeoJSON points, set a `_score` value in your return table and `declutter_below` on the layer. + Shapefiles/GeoJSON **must** be in WGS84 projection, i.e. pure latitude/longitude. (Use ogr2ogr to reproject them if your source material is in a different projection.) They will be clipped to the bounds of the first .pbf that you import, unless you specify otherwise with a `bounding_box` setting in your JSON file. ### Lua spatial queries diff --git a/include/declutter.h b/include/declutter.h new file mode 100644 index 00000000..7b4201fe --- /dev/null +++ b/include/declutter.h @@ -0,0 +1,56 @@ +/*! \file */ +#ifndef _DECLUTTER_H +#define _DECLUTTER_H + +#include +#include +#include "coordinates.h" +#include "output_object.h" + +class TileDataSource; +struct LayerDef; + +/// A point feature held back from the tile index until its minimum zoom has been decided +struct DeclutterEntry { + OutputObject oo; + LatpLon point; + uint64_t id; + int32_t score; + bool fromShapefile; +}; + +/** + \brief Thins out point features so only the most important ones appear at low zooms. + + Layers with declutter_below set don't index their point features as they're read: the + features are parked here, together with the score their Lua profile gave them with + Score(). Once everything has been read, apply() works up through the zoom levels, + placing the highest-scoring features first and holding back any that fall too close to + one already placed, then hands them all to the tile index. + + Ranking globally rather than per z6 tile means features either side of a z6 boundary + still compete with each other. +*/ +class Declutter { + +public: + /// Note which layers are decluttered - call before any add() + void configure(const std::vector& layers); + + bool inUse() const { return anyLayers; } + bool isDecluttered(uint_least8_t layer) const { return anyLayers && layerEnabled[layer]; } + + /// Park a point feature (thread-safe) + void add(const OutputObject& oo, LatpLon point, uint64_t id, int32_t score, bool fromShapefile); + + /// Assign minimum zooms, then write everything to the tile index + void apply(const std::vector& layers, TileDataSource& osmSource, TileDataSource& shpSource); + +private: + bool anyLayers = false; + std::vector layerEnabled; + std::vector> entries; // one per layer + std::vector mutexes; // | +}; + +#endif //_DECLUTTER_H diff --git a/include/geojson_processor.h b/include/geojson_processor.h index 51d45d01..442dc2f3 100644 --- a/include/geojson_processor.h +++ b/include/geojson_processor.h @@ -43,7 +43,7 @@ class GeoJSONProcessor { template std::vector pointsFromGeoJSONArray(const rapidjson::GenericArray &arr); - AttributeIndex readProperties(const rapidjson::Value &pr, bool &hasName, std::string &name, LayerDef &layer, unsigned &minzoom); + AttributeIndex readProperties(const rapidjson::Value &pr, bool &hasName, std::string &name, LayerDef &layer, unsigned &minzoom, int32_t &score); }; #endif //_GEOJSON_PROCESSOR_H diff --git a/include/osm_lua_processing.h b/include/osm_lua_processing.h index b8d536c1..912652d8 100644 --- a/include/osm_lua_processing.h +++ b/include/osm_lua_processing.h @@ -11,6 +11,7 @@ #include "osm_store.h" #include "shared_data.h" #include "output_object.h" +#include "declutter.h" #include "shp_mem_tiles.h" #include "osm_mem_tiles.h" #include "helpers.h" @@ -57,6 +58,7 @@ class OsmLuaProcessing { const class ShpMemTiles &shpMemTiles, class OsmMemTiles &osmMemTiles, AttributeStore &attributeStore, + class Declutter &declutter, bool materializeGeometries, bool isFirst ); @@ -204,6 +206,7 @@ class OsmLuaProcessing { void AttributeInteger(const std::string &key, const int val, const char minzoom); void MinZoom(const double z); void ZOrder(const double z); + void Score(const double score); // Relation scan support @@ -265,6 +268,7 @@ class OsmLuaProcessing { relationAccepted = false; relationList.clear(); relationSubscript = -1; + declutterOutputs.clear(); lastStoredGeometryId = 0; isWay = false; isRelation = false; @@ -273,6 +277,8 @@ class OsmLuaProcessing { void removeAttributeIfNeeded(const std::string& key); + void noteIfDecluttered(LatpLon point); + const inline Point getPoint() { return Point(lon/10000000.0,latp/10000000.0); } @@ -291,6 +297,7 @@ class OsmLuaProcessing { const class ShpMemTiles &shpMemTiles; class OsmMemTiles &osmMemTiles; AttributeStore &attributeStore; // key/value store + class Declutter &declutter; // point features held back for decluttering int64_t originalOsmID; ///< Original OSM object ID bool isWay, isRelation, isClosed; ///< Way, node, relation? @@ -321,6 +328,11 @@ class OsmLuaProcessing { class LayerDefinition &layers; std::vector> outputs; // All output objects that have been created + + // Point features in decluttered layers: the exact location we'd have indexed them at, + // plus the score from Score(). Sparse - most objects add nothing here. + struct PendingDeclutter { uint32_t outputIndex; LatpLon point; int32_t score; }; + std::vector declutterOutputs; std::vector outputKeys; const PbfReader::Relation* currentRelation; const boost::container::flat_map* currentPostScanTags; // for postScan only diff --git a/include/shared_data.h b/include/shared_data.h index d367a256..5cfe13af 100644 --- a/include/shared_data.h +++ b/include/shared_data.h @@ -39,6 +39,10 @@ struct LayerDef { std::string indexName; std::map attributeMap; // string 0, number 1, bool 2 bool writeTo; + // Decluttering: set after addLayer(), so these live at the end of the struct + uint declutterBelow = 0; // zoom below which point features are thinned out by Score() + double declutterDistance = 40; // how far apart (in 256px screen pixels) they should be kept + double declutterThreshold = 0; // score needed at the layer's minzoom (halves at each zoom) const bool useColumn(std::string &col) { return allSourceColumns || (std::find(sourceColumns.begin(), sourceColumns.end(), col) != sourceColumns.end()); diff --git a/include/shp_mem_tiles.h b/include/shp_mem_tiles.h index 80ae4496..a5cf4f3f 100644 --- a/include/shp_mem_tiles.h +++ b/include/shp_mem_tiles.h @@ -3,13 +3,14 @@ #define _SHP_MEM_TILES #include "tile_data.h" +#include "declutter.h" extern bool verbose; class ShpMemTiles : public TileDataSource { public: - ShpMemTiles(size_t threadNum, uint indexZoom); + ShpMemTiles(size_t threadNum, uint indexZoom, class Declutter& declutter); std::string name() const override { return "shp"; } @@ -25,6 +26,7 @@ class ShpMemTiles : public TileDataSource bool hasName, const std::string& name, uint minzoom, + int32_t score, AttributeIndex attrIdx ); @@ -62,6 +64,7 @@ class ShpMemTiles : public TileDataSource } private: + class Declutter& declutter; std::vector indexedGeometries; // prepared boost::geometry objects (from shapefiles) std::map indexedGeometryNames; // | optional names for each one std::map indices; // Spatial indices, boost::geometry::index objects for shapefile indices diff --git a/include/shp_processor.h b/include/shp_processor.h index 386ab93b..14aeeee1 100644 --- a/include/shp_processor.h +++ b/include/shp_processor.h @@ -41,11 +41,11 @@ class ShpProcessor { AttributeIndex readShapefileAttributes(DBFHandle dbf, int recordNum, std::unordered_map &columnMap, std::unordered_map &columnTypeMap, - LayerDef &layer, uint &minzoom); + LayerDef &layer, uint &minzoom, int32_t &score); // Process an individual shapefile record void processShapeGeometry(SHPObject* shape, AttributeIndex attrIdx, - const LayerDef &layer, uint layerNum, bool hasName, const std::string &name); + const LayerDef &layer, uint layerNum, bool hasName, const std::string &name, int32_t score); }; #endif //_SHP_PROCESSOR_H diff --git a/resources/config-schema.json b/resources/config-schema.json index 34a9d1c4..ca8841c7 100644 --- a/resources/config-schema.json +++ b/resources/config-schema.json @@ -59,6 +59,9 @@ "simplify_ratio": { "type": "number" }, "filter_below": { "type": "integer", "minimum": 0 }, "filter_area": { "type": "number" }, + "declutter_below": { "type": "integer", "minimum": 0 }, + "declutter_distance": { "type": "number" }, + "declutter_threshold": { "type": "number" }, "feature_limit": { "type": "integer", "minimum": 0 }, "feature_limit_below": { "type": "integer", "minimum": 0 }, "combine_points": { "type": "boolean" }, diff --git a/src/declutter.cpp b/src/declutter.cpp new file mode 100644 index 00000000..d41cd5f3 --- /dev/null +++ b/src/declutter.cpp @@ -0,0 +1,110 @@ +#include "declutter.h" + +#include +#include +#include + +#include "shared_data.h" +#include "tile_data.h" + +namespace bgi = boost::geometry::index; + +// Features are ranked in world units (0-1 across the map in each direction), so a +// separation of n screen pixels at zoom z is simply n/(256< WorldPoint; + +void Declutter::configure(const std::vector& layers) { + layerEnabled.assign(layers.size(), false); + for (size_t i = 0; i < layers.size(); i++) { + if (layers[i].declutterBelow == 0) continue; + layerEnabled[i] = true; + anyLayers = true; + } + if (!anyLayers) return; + entries.resize(layers.size()); + mutexes = std::vector(layers.size()); +} + +void Declutter::add(const OutputObject& oo, LatpLon point, uint64_t id, int32_t score, bool fromShapefile) { + std::lock_guard lock(mutexes[oo.layer]); + entries[oo.layer].push_back({ oo, point, id, score, fromShapefile }); +} + +// Work up through the zoom levels, giving each feature the lowest zoom at which it both +// clears the score threshold and isn't crowded out by a feature already placed. The +// threshold halves at each zoom, and the separation is constant in screen terms, so each +// zoom admits both lower-scoring and more tightly packed features than the one before. +static void assignMinZooms(const LayerDef& layer, std::vector& list) { + // Highest score first; ties broken on id and position so the ranking doesn't depend on + // the order features happened to be read in + std::sort(list.begin(), list.end(), [](const DeclutterEntry& a, const DeclutterEntry& b) { + if (a.score != b.score) return a.score > b.score; + if (a.id != b.id) return a.id < b.id; + if (a.point.latp != b.point.latp) return a.point.latp < b.point.latp; + return a.point.lon < b.point.lon; + }); + + bgi::rtree> placed; + std::vector done(list.size(), false); + size_t remaining = list.size(); + + double threshold = layer.declutterThreshold; + for (uint z = layer.minzoom; z < layer.declutterBelow && remaining > 0; z++, threshold /= 2) { + const double separation = layer.declutterDistance / (256.0 * (1u << z)); + + for (size_t i = 0; i < list.size(); i++) { + DeclutterEntry& e = list[i]; + if (e.score < threshold) break; // sorted by descending score, so nothing further qualifies + if (done[i] || e.oo.minZoom > z) continue; + + // Tweak for the Cheltenham case (two significant cities next to each other) - don't let it be repeatedly pushed out by Gloucester + double allowed = separation; + if (threshold > 0 && e.score > threshold * 5) allowed = separation / (e.score / threshold / 3.0); + + const double x = lon2tilexf(e.point.lon / 10000000.0, 0); + const double y = latp2tileyf(e.point.latp / 10000000.0, 0); + const WorldPoint p(x, y); + + // Anything within `allowed` of p is inside this box, so we only have to measure + // the handful of features the box picks up + const boost::geometry::model::box around( + WorldPoint(x - allowed, y - allowed), WorldPoint(x + allowed, y + allowed)); + bool crowded = false; + for (auto it = placed.qbegin(bgi::intersects(around)); it != placed.qend() && !crowded; ++it) + crowded = boost::geometry::distance(p, *it) < allowed; + if (crowded) continue; + + e.oo.setMinZoom(z); + done[i] = true; + remaining--; + placed.insert(p); + } + } + + // Anything that never won a place appears from declutter_below upwards + for (size_t i = 0; i < list.size(); i++) + if (!done[i] && list[i].oo.minZoom < layer.declutterBelow) + list[i].oo.setMinZoom(layer.declutterBelow); +} + +void Declutter::apply(const std::vector& layers, TileDataSource& osmSource, TileDataSource& shpSource) { + if (!anyLayers) return; + + for (size_t i = 0; i < entries.size(); i++) { + if (entries[i].empty()) continue; + std::cout << "Decluttering " << entries[i].size() << " features in layer " << layers[i].name << ":" << std::flush; + assignMinZooms(layers[i], entries[i]); + + std::vector perZoom(16, 0); // minZoom is a 4-bit field + for (const auto& e : entries[i]) perZoom[e.oo.minZoom]++; + for (size_t z = 0; z < perZoom.size(); z++) + if (perZoom[z] > 0) std::cout << " z" << z << ":" << perZoom[z]; + std::cout << std::endl; + + for (const auto& e : entries[i]) { + TileDataSource& source = e.fromShapefile ? shpSource : osmSource; + source.addObjectToSmallIndex(latpLon2index(e.point, source.getIndexZoom()), e.oo, e.id); + } + std::vector().swap(entries[i]); + } +} diff --git a/src/geojson_processor.cpp b/src/geojson_processor.cpp index 131458f7..fa493ac0 100644 --- a/src/geojson_processor.cpp +++ b/src/geojson_processor.cpp @@ -112,7 +112,8 @@ void GeoJSONProcessor::processFeature(rapidjson::GenericObject feature, std::string name; const rapidjson::Value &pr = feature["properties"]; unsigned minzoom = layer.minzoom; - AttributeIndex attrIdx = readProperties(pr, hasName, name, layer, minzoom); + int32_t score = 0; + AttributeIndex attrIdx = readProperties(pr, hasName, name, layer, minzoom, score); // Parse geometry auto geometry = feature["geometry"].GetObject(); @@ -130,7 +131,7 @@ void GeoJSONProcessor::processFeature(rapidjson::GenericObject feature, // coordinates is [x,y] Point p( coords[0].GetDouble(), lat2latp(coords[1].GetDouble()) ); if (geom::within(p, clippingBox)) { - shpMemTiles.StoreGeometry(layerNum, layer.name, POINT_, p, layer.indexed, hasName, name, minzoom, attrIdx); + shpMemTiles.StoreGeometry(layerNum, layer.name, POINT_, p, layer.indexed, hasName, name, minzoom, score, attrIdx); } } else if (geomType=="LineString") { @@ -140,7 +141,7 @@ void GeoJSONProcessor::processFeature(rapidjson::GenericObject feature, MultiLinestring out; geom::intersection(ls, clippingBox, out); if (!geom::is_empty(out)) { - shpMemTiles.StoreGeometry(layerNum, layer.name, MULTILINESTRING_, out, layer.indexed, hasName, name, minzoom, attrIdx); + shpMemTiles.StoreGeometry(layerNum, layer.name, MULTILINESTRING_, out, layer.indexed, hasName, name, minzoom, score, attrIdx); } } else if (geomType=="Polygon") { @@ -151,7 +152,7 @@ void GeoJSONProcessor::processFeature(rapidjson::GenericObject feature, MultiPolygon out; geom::intersection(polygon, clippingBox, out); if (!geom::is_empty(out)) { - shpMemTiles.StoreGeometry(layerNum, layer.name, POLYGON_, out, layer.indexed, hasName, name, minzoom, attrIdx); + shpMemTiles.StoreGeometry(layerNum, layer.name, POLYGON_, out, layer.indexed, hasName, name, minzoom, score, attrIdx); } } else if (geomType=="MultiPoint") { @@ -159,7 +160,7 @@ void GeoJSONProcessor::processFeature(rapidjson::GenericObject feature, for (auto &pt : coords) { Point p( pt[0].GetDouble(), lat2latp(pt[1].GetDouble()) ); if (geom::within(p, clippingBox)) { - shpMemTiles.StoreGeometry(layerNum, layer.name, POINT_, p, layer.indexed, hasName, name, minzoom, attrIdx); + shpMemTiles.StoreGeometry(layerNum, layer.name, POINT_, p, layer.indexed, hasName, name, minzoom, score, attrIdx); } } @@ -174,7 +175,7 @@ void GeoJSONProcessor::processFeature(rapidjson::GenericObject feature, MultiLinestring out; geom::intersection(mls, clippingBox, out); if (!geom::is_empty(out)) { - shpMemTiles.StoreGeometry(layerNum, layer.name, MULTILINESTRING_, out, layer.indexed, hasName, name, minzoom, attrIdx); + shpMemTiles.StoreGeometry(layerNum, layer.name, MULTILINESTRING_, out, layer.indexed, hasName, name, minzoom, score, attrIdx); } } else if (geomType=="MultiPolygon") { @@ -187,7 +188,7 @@ void GeoJSONProcessor::processFeature(rapidjson::GenericObject feature, MultiPolygon out; geom::intersection(mp, clippingBox, out); if (!geom::is_empty(out)) { - shpMemTiles.StoreGeometry(layerNum, layer.name, POLYGON_, out, layer.indexed, hasName, name, minzoom, attrIdx); + shpMemTiles.StoreGeometry(layerNum, layer.name, POLYGON_, out, layer.indexed, hasName, name, minzoom, score, attrIdx); } } } @@ -215,7 +216,7 @@ std::vector GeoJSONProcessor::pointsFromGeoJSONArray(const rapidjson::Gen } // Read properties and generate an AttributeIndex -AttributeIndex GeoJSONProcessor::readProperties(const rapidjson::Value &pr, bool &hasName, std::string &name, LayerDef &layer, unsigned &minzoom) { +AttributeIndex GeoJSONProcessor::readProperties(const rapidjson::Value &pr, bool &hasName, std::string &name, LayerDef &layer, unsigned &minzoom, int32_t &score) { std::lock_guard lock(attributeMutex); AttributeStore& attributeStore = osmLuaProcessing.getAttributeStore(); AttributeSet attributes; @@ -259,6 +260,7 @@ AttributeIndex GeoJSONProcessor::readProperties(const rapidjson::Value &pr, bool layer.attributeMap[key] = 0; } else if (val.isType()) { if (key=="_minzoom") { minzoom=val; continue; } + if (key=="_score") { score=val; continue; } attributeStore.addAttribute(attributes, key, (int)val, 0); layer.attributeMap[key] = 1; } else if (val.isType()) { diff --git a/src/osm_lua_processing.cpp b/src/osm_lua_processing.cpp index 62303238..b95b56d2 100644 --- a/src/osm_lua_processing.cpp +++ b/src/osm_lua_processing.cpp @@ -194,6 +194,7 @@ void rawLayer(const std::string& layerName, bool area) { return osmLuaProcessing void rawLayerAsCentroid(const std::string &layerName, kaguya::VariadicArgType nodeSources) { return osmLuaProcessing->LayerAsCentroid(layerName, nodeSources); } void rawMinZoom(const double z) { return osmLuaProcessing->MinZoom(z); } void rawZOrder(const double z) { return osmLuaProcessing->ZOrder(z); } +void rawScore(const double score) { return osmLuaProcessing->Score(score); } OsmLuaProcessing::OptionalRelation rawNextRelation() { return osmLuaProcessing->NextRelation(); } void rawRestartRelations() { return osmLuaProcessing->RestartRelations(); } std::string rawFindInRelation(const std::string& key) { return osmLuaProcessing->FindInRelation(key); } @@ -230,12 +231,14 @@ OsmLuaProcessing::OsmLuaProcessing( const class ShpMemTiles &shpMemTiles, class OsmMemTiles &osmMemTiles, AttributeStore &attributeStore, + class Declutter &declutter, bool materializeGeometries, bool isFirst) : osmStore(osmStore), shpMemTiles(shpMemTiles), osmMemTiles(osmMemTiles), attributeStore(attributeStore), + declutter(declutter), config(configIn), currentTags(NULL), layers(layers), @@ -292,6 +295,7 @@ OsmLuaProcessing::OsmLuaProcessing( luaState["MinZoom"] = &rawMinZoom; luaState["ZOrder"] = &rawZOrder; + luaState["Score"] = &rawScore; luaState["Accept"] = &rawAccept; luaState["NextRelation"] = &rawNextRelation; luaState["RestartRelations"] = &rawRestartRelations; @@ -657,6 +661,7 @@ void OsmLuaProcessing::Layer(const string &layerName, bool area) { id = osmMemTiles.storePoint(p); OutputObject oo(geomType, layers.layerMap[layerName], id, 0, layerMinZoom); outputs.push_back(std::make_pair(std::move(oo), attributes)); + noteIfDecluttered(LatpLon { latp, lon }); return; } else if (geomType==POLYGON_) { @@ -860,6 +865,14 @@ void OsmLuaProcessing::LayerAsCentroid(const string &layerName, kaguya::Variadic } OutputObject oo(POINT_, layers.layerMap[layerName], id, 0, layerMinZoom); outputs.push_back(std::make_pair(std::move(oo), attributes)); + noteIfDecluttered(LatpLon { (int32_t)geomp.y(), (int32_t)geomp.x() }); +} + +// If the layer we've just written a point to is decluttered, remember where the point is: +// it'll be held back from the tile index until all features have been read and ranked. +void OsmLuaProcessing::noteIfDecluttered(LatpLon point) { + if (!declutter.isDecluttered(outputs.back().first.layer)) return; + declutterOutputs.push_back({ (uint32_t)(outputs.size() - 1), point, 0 }); } Point OsmLuaProcessing::calculateCentroid(CentroidAlgorithm algorithm) { @@ -1005,6 +1018,16 @@ void OsmLuaProcessing::MinZoom(const double z) { outputs.back().first.setMinZoom(z); } +// Set the score used to rank this feature against its neighbours when decluttering +void OsmLuaProcessing::Score(const double score) { + if (outputs.size()==0) { ProcessingError("Can't set score if no Layer set"); return; } + if (declutterOutputs.empty() || declutterOutputs.back().outputIndex != outputs.size()-1) { + ProcessingError("Score() only applies to point features in a layer with declutter_below set"); + return; + } + declutterOutputs.back().score = OutputObject::finite_cast(score); +} + // Set z_order void OsmLuaProcessing::ZOrder(const double z) { if (outputs.size()==0) { ProcessingError("Can't set z_order if no Layer set"); return; } @@ -1179,7 +1202,8 @@ bool OsmLuaProcessing::setWay(WayID wayId, LatpLonVec const &llVec, const TagMap } if (!this->empty()) { - osmMemTiles.addGeometryToIndex(linestringCached(), finalizeOutputs(), originalOsmID); + std::vector objects = finalizeOutputs(); + if (!objects.empty()) osmMemTiles.addGeometryToIndex(linestringCached(), objects, originalOsmID); return wayEmitted; } @@ -1227,11 +1251,12 @@ void OsmLuaProcessing::setRelation( if (this->empty()) return; try { + std::vector objects = finalizeOutputs(); + if (objects.empty()) return; if (isClosed) { - std::vector objects = finalizeOutputs(); osmMemTiles.addGeometryToIndex(multiPolygonCached(), objects, originalOsmID); } else { - osmMemTiles.addGeometryToIndex(multiLinestringCached(), finalizeOutputs(), originalOsmID); + osmMemTiles.addGeometryToIndex(multiLinestringCached(), objects, originalOsmID); } } catch(std::out_of_range &err) { cout << "In relation " << originalOsmID << ": " << err.what() << endl; @@ -1260,9 +1285,17 @@ SignificantTags OsmLuaProcessing::GetSignificantWayKeys() { std::vector OsmLuaProcessing::finalizeOutputs() { std::vector list; list.reserve(this->outputs.size()); - for (auto jt = this->outputs.begin(); jt != this->outputs.end(); ++jt) { - jt->first.setAttributeSet(attributeStore.add(jt->second)); - list.push_back(jt->first); + size_t nextDeclutter = 0; + for (size_t i = 0; i < this->outputs.size(); i++) { + auto& output = this->outputs[i]; + output.first.setAttributeSet(attributeStore.add(output.second)); + if (nextDeclutter < declutterOutputs.size() && declutterOutputs[nextDeclutter].outputIndex == i) { + // Park it: Declutter::apply() will set its minimum zoom and index it later + const auto& pending = declutterOutputs[nextDeclutter++]; + declutter.add(output.first, pending.point, originalOsmID, pending.score, false); + continue; + } + list.push_back(output.first); } return list; } diff --git a/src/shared_data.cpp b/src/shared_data.cpp index 60b6979e..d3843562 100644 --- a/src/shared_data.cpp +++ b/src/shared_data.cpp @@ -343,12 +343,25 @@ void Config::readConfig(rapidjson::Document &jsonConfig, bool &hasClippingBox, B } string indexName = it->value.HasMember("index_column") ? it->value["index_column"].GetString() : ""; - layers.addLayer(layerName, minZoom, maxZoom, + uint layerNum = layers.addLayer(layerName, minZoom, maxZoom, simplifyBelow, simplifyLevel, simplifyLength, simplifyRatio, simplifyAlgo, filterBelow, filterArea, sortZOrderAscending, featureLimit, featureLimitBelow, combinePoints, combineLinesBelow, combinePolyBelow, source, sourceColumns, allSourceColumns, indexed, indexName, writeTo); + // Decluttering (thinning out point features by the score their profile gives them) + if (it->value.HasMember("declutter_below")) { + LayerDef &layerDef = layers.layers[layerNum]; + int declutterBelow = it->value["declutter_below"].GetInt(); + if (declutterBelow > 15) { // minZoom is a 4-bit field + cerr << "declutter_below in layer " << layerName << " capped at z15" << endl; + declutterBelow = 15; + } + layerDef.declutterBelow = declutterBelow; + if (it->value.HasMember("declutter_distance" )) layerDef.declutterDistance = it->value["declutter_distance" ].GetDouble(); + if (it->value.HasMember("declutter_threshold")) layerDef.declutterThreshold = it->value["declutter_threshold"].GetDouble(); + } + cout << "Layer " << layerName << " (z" << minZoom << "-" << maxZoom << ")"; if (it->value.HasMember("write_to")) { cout << " -> " << it->value["write_to"].GetString(); } cout << endl; diff --git a/src/shp_mem_tiles.cpp b/src/shp_mem_tiles.cpp index b417423e..977a37da 100644 --- a/src/shp_mem_tiles.cpp +++ b/src/shp_mem_tiles.cpp @@ -7,8 +7,9 @@ using namespace std; namespace geom = boost::geometry; extern bool verbose; -ShpMemTiles::ShpMemTiles(size_t threadNum, uint indexZoom) +ShpMemTiles::ShpMemTiles(size_t threadNum, uint indexZoom, class Declutter& declutter) : TileDataSource(threadNum, indexZoom, false), + declutter(declutter), spatialIndexZoom(15) { } @@ -135,6 +136,7 @@ void ShpMemTiles::StoreGeometry( bool hasName, const std::string& name, uint minzoom, + int32_t score, AttributeIndex attrIdx ) { @@ -151,9 +153,14 @@ void ShpMemTiles::StoreGeometry( Point sp(p->x()*10000000.0, p->y()*10000000.0); NodeID oid = storePoint(sp); oo = std::make_shared(geomType, layerNum, oid, attrIdx, minzoom); - tilex = lon2tilex(p->x(), indexZoom); - tiley = latp2tiley(p->y(), indexZoom); - addObjectToSmallIndex(TileCoordinates(tilex, tiley), *oo, 0); + if (declutter.isDecluttered(layerNum)) { + // held back until all features have been read and ranked + declutter.add(*oo, LatpLon { (int32_t)sp.y(), (int32_t)sp.x() }, 0, score, true); + } else { + tilex = lon2tilex(p->x(), indexZoom); + tiley = latp2tiley(p->y(), indexZoom); + addObjectToSmallIndex(TileCoordinates(tilex, tiley), *oo, 0); + } } else { return; } } break; diff --git a/src/shp_processor.cpp b/src/shp_processor.cpp index 469932c8..4cae1179 100644 --- a/src/shp_processor.cpp +++ b/src/shp_processor.cpp @@ -58,7 +58,7 @@ void ShpProcessor::fillPointArrayFromShapefile(vector *points, SHPObject AttributeIndex ShpProcessor::readShapefileAttributes( DBFHandle dbf, int recordNum, unordered_map &columnMap, unordered_map &columnTypeMap, - LayerDef &layer, uint &minzoom) { + LayerDef &layer, uint &minzoom, int32_t &score) { std::lock_guard lock(attributeMutex); AttributeStore& attributeStore = osmLuaProcessing.getAttributeStore(); @@ -89,6 +89,7 @@ AttributeIndex ShpProcessor::readShapefileAttributes( layer.attributeMap[key] = 0; } else if (val.isType()) { if (key=="_minzoom") { minzoom=val; continue; } + if (key=="_score") { score=val; continue; } attributeStore.addAttribute(attributes, key, (int)val, 0); layer.attributeMap[key] = 1; } else if (val.isType()) { @@ -188,9 +189,10 @@ void ShpProcessor::read(class LayerDef &layer, uint layerNum) std::lock_guard lock(attributeMutex); name=DBFReadStringAttribute(dbf.get(), i, indexField); hasName = true; } - AttributeIndex attrIdx = readShapefileAttributes(dbf.get(), i, columnMap, columnTypeMap, layer, layer.minzoom); + int32_t score = 0; + AttributeIndex attrIdx = readShapefileAttributes(dbf.get(), i, columnMap, columnTypeMap, layer, layer.minzoom, score); // process geometry - processShapeGeometry(shape.get(), attrIdx, layer, layerNum, hasName, name); + processShapeGeometry(shape.get(), attrIdx, layer, layerNum, hasName, name, score); } catch (...) { std::lock_guard lock(errorMutex); if (!error) @@ -204,7 +206,7 @@ void ShpProcessor::read(class LayerDef &layer, uint layerNum) } void ShpProcessor::processShapeGeometry(SHPObject* shape, AttributeIndex attrIdx, - const LayerDef &layer, uint layerNum, bool hasName, const string &name) { + const LayerDef &layer, uint layerNum, bool hasName, const string &name, int32_t score) { int shapeType = shape->nSHPType; // 1=point, 3=polyline, 5=(multi)polygon [8=multipoint, 11+=3D] int minzoom = layer.minzoom; @@ -212,7 +214,7 @@ void ShpProcessor::processShapeGeometry(SHPObject* shape, AttributeIndex attrIdx // Points Point p( shape->padfX[0], lat2latp(shape->padfY[0]) ); if (geom::within(p, clippingBox)) { - shpMemTiles.StoreGeometry(layerNum, layer.name, POINT_, p, layer.indexed, hasName, name, minzoom, attrIdx); + shpMemTiles.StoreGeometry(layerNum, layer.name, POINT_, p, layer.indexed, hasName, name, minzoom, score, attrIdx); } } else if (shapeType==8 || shapeType==18 || shapeType==28) { @@ -220,7 +222,7 @@ void ShpProcessor::processShapeGeometry(SHPObject* shape, AttributeIndex attrIdx for (uint i=0; inVertices; i++) { Point p( shape->padfX[i], lat2latp(shape->padfY[i]) ); if (geom::within(p, clippingBox)) { - shpMemTiles.StoreGeometry(layerNum, layer.name, POINT_, p, layer.indexed, hasName, name, minzoom, attrIdx); + shpMemTiles.StoreGeometry(layerNum, layer.name, POINT_, p, layer.indexed, hasName, name, minzoom, score, attrIdx); } } @@ -236,7 +238,7 @@ void ShpProcessor::processShapeGeometry(SHPObject* shape, AttributeIndex attrIdx MultiLinestring out; geom::intersection(ls, clippingBox, out); for (MultiLinestring::const_iterator it = out.begin(); it != out.end(); ++it) { - shpMemTiles.StoreGeometry(layerNum, layer.name, LINESTRING_, *it, layer.indexed, hasName, name, minzoom, attrIdx); + shpMemTiles.StoreGeometry(layerNum, layer.name, LINESTRING_, *it, layer.indexed, hasName, name, minzoom, score, attrIdx); } } @@ -295,7 +297,7 @@ void ShpProcessor::processShapeGeometry(SHPObject* shape, AttributeIndex attrIdx MultiPolygon out; geom::intersection(multi, clippingBox, out); if (boost::size(out)>0) { - shpMemTiles.StoreGeometry(layerNum, layer.name, POLYGON_, out, layer.indexed, hasName, name, minzoom, attrIdx); + shpMemTiles.StoreGeometry(layerNum, layer.name, POLYGON_, out, layer.indexed, hasName, name, minzoom, score, attrIdx); } } else { diff --git a/src/tilemaker.cpp b/src/tilemaker.cpp index 79a3ed8d..f4891425 100644 --- a/src/tilemaker.cpp +++ b/src/tilemaker.cpp @@ -261,14 +261,17 @@ int main(const int argc, const char* argv[]) { class LayerDefinition layers(config.layers); + class Declutter declutter; + declutter.configure(layers.layers); + const unsigned int indexZoom = std::min(config.baseZoom, 14u); class OsmMemTiles osmMemTiles(options.threadNum, indexZoom, config.includeID, *nodeStore, *wayStore); - class ShpMemTiles shpMemTiles(options.threadNum, indexZoom); + class ShpMemTiles shpMemTiles(options.threadNum, indexZoom, declutter); osmMemTiles.open(); shpMemTiles.open(); OsmLuaProcessing osmLuaProcessing(osmStore, config, layers, options.luaFile, - shpMemTiles, osmMemTiles, attributeStore, options.osm.materializeGeometries, true); + shpMemTiles, osmMemTiles, attributeStore, declutter, options.osm.materializeGeometries, true); // ---- Load external sources (shp/geojson) @@ -336,7 +339,7 @@ int main(const int argc, const char* argv[]) { [&]() { thread_local std::pair> osmLuaProcessing; if (osmLuaProcessing.first != inputFile) { - osmLuaProcessing = std::make_pair(inputFile, std::make_shared(osmStore, config, layers, options.luaFile, shpMemTiles, osmMemTiles, attributeStore, options.osm.materializeGeometries, false)); + osmLuaProcessing = std::make_pair(inputFile, std::make_shared(osmStore, config, layers, options.luaFile, shpMemTiles, osmMemTiles, attributeStore, declutter, options.osm.materializeGeometries, false)); } return osmLuaProcessing.second; }, @@ -378,6 +381,9 @@ int main(const int argc, const char* argv[]) { // Loop through tiles std::atomic tilesWritten(0), lastTilesWritten(0); + // Rank any decluttered point features, and index them at their allotted zoom levels + declutter.apply(layers.layers, osmMemTiles, shpMemTiles); + for (auto source : sources) { source->finalize(options.threadNum); }