[Git][debian-gis-team/tilemaker][upstream] New upstream version 2.4.0

ǝɹʇʇɐʃǝ◖ xıʃǝɟ (@pantierra) gitlab at salsa.debian.org
Sat Jul 1 09:56:05 BST 2023



ǝɹʇʇɐʃǝ◖ xıʃǝɟ pushed to branch upstream at Debian GIS Project / tilemaker


Commits:
1ca89b0c by Felix Delattre at 2023-07-01T08:44:02+00:00
New upstream version 2.4.0
- - - - -


20 changed files:

- CHANGELOG.md
- README.md
- docs/CONFIGURATION.md
- docs/INSTALL.md
- include/coordinates.h
- include/geom.h
- include/osm_lua_processing.h
- include/output_object.h
- include/shared_data.h
- include/shp_mem_tiles.h
- include/tile_data.h
- resources/process-openmaptiles.lua
- src/geom.cpp
- src/osm_lua_processing.cpp
- src/output_object.cpp
- src/read_shp.cpp
- src/shared_data.cpp
- src/shp_mem_tiles.cpp
- src/tile_data.cpp
- src/tilemaker.cpp


Changes:

=====================================
CHANGELOG.md
=====================================
@@ -1,6 +1,24 @@
 # Changelog
 
-# [2.3.0] - 2023-03-08
+
+## [2.4.0] - 2023-03-28
+
+### Added
+- Option to reverse object sort order (@Nakaner, @systemed)
+- Compile-time option to use floats for ZOrder (@Nakaner, @systemed)
+- Advisory note if user tries to generate tiles at z16+ (@systemed)
+
+### Changed
+- Faster tile clipping (@systemed based on code by @mourner)
+- Use rtree to index large polygons (@systemed, @kleunen)
+
+### Fixed
+- Update use of access in OpenMapTiles-compatible schema (@dschep)
+- Align path/track transportation classes with OpenMapTiles (@dschep)
+- Add missing paved/unpaved values as per OpenMapTiles (@dschep)
+
+
+## [2.3.0] - 2023-03-08
 
 ### Added
 - Send project name to init_function (@systemed)


=====================================
README.md
=====================================
@@ -90,4 +90,4 @@ Formatting: braces and indents as shown, hard tabs (4sp). (Yes, I know.) Please
 
 Tilemaker is maintained by Richard Fairhurst and supported by [many contributors](https://github.com/systemed/tilemaker/graphs/contributors).
 
-Copyright tilemaker contributors, 2015-2021. The tilemaker code is licensed as FTWPL; you may do anything you like with this code and there is no warranty. The included sqlite_modern_cpp (Amin Roosta) is MIT; [kaguya](https://github.com/satoren/kaguya) is licensed under the Boost Software Licence.
+Copyright tilemaker contributors, 2015-2023. The tilemaker code is licensed as FTWPL; you may do anything you like with this code and there is no warranty. The included sqlite_modern_cpp (Amin Roosta) is MIT; [kaguya](https://github.com/satoren/kaguya) is licensed under the Boost Software Licence.


=====================================
docs/CONFIGURATION.md
=====================================
@@ -73,6 +73,7 @@ You can add optional parameters to layers:
 * `filter_below` - filter areas by minimum size below this zoom level
 * `filter_area` - minimum size (in square degrees of longitude) for the zoom level `filter_below-1`
 * `combine_polygons_below` - merge adjacent polygons with the same attributes below this zoom level
+* `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)
 
 Use these options to combine different layer specs within one outputted layer. For example:
 
@@ -133,7 +134,7 @@ To do that, you use these methods:
 * `node:Attribute(key,value)` or `node:Attribute(key,value)`: add an attribute to the most recently written layer.
 * `node:AttributeNumeric(key,value)`, `node:AttributeBoolean(key,value)` (and `way:`...): for numeric/boolean columns.
 * `node:Id()` or `way:Id()`: get the OSM ID of the current object.
-* `node:ZOrder(number)` or `way:ZOrder(number)`: Set a numeric value (default 0, 1-byte unsigned integer) 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` or `combine_polygons_below` is used.
+* `node:ZOrder(number)` or `way:ZOrder(number)`: Set a numeric value (default 0, 1-byte signed integer) 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` or `combine_polygons_below` is used.
 * `node:MinZoom(zoom)` or `way: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).
 * `way:Length()` and `way:Area()`: return the length (metres)/area (square metres) of the current object. Requires recent Boost.
 * `way: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).


=====================================
docs/INSTALL.md
=====================================
@@ -63,3 +63,13 @@ The docker container can be run like this:
     docker run -v /Users/Local/Downloads/:/srv -i -t --rm tilemaker /srv/germany-latest.osm.pbf --output=/srv/germany.mbtiles
 
 Keep in mind to map the volume your .osm.pbf files are in to a path within your docker container, as seen in the example above. 
+
+### Compile-time options
+
+tilemaker has two compile-time options that increase memory usage but may be useful in certain circumstances. You can include them when building like this:
+
+    make "CONFIG=-DFLOAT_Z_ORDER"
+
+FLOAT_Z_ORDER allows you to use a full range of ZOrder values in your Lua script, rather than being restricted to single-byte integer (-127 to 127).
+
+FAT_TILE_INDEX allows you to generate vector tiles at zoom level 17 or greater. You almost certainly don't need to do this. Vector tiles are usually generated up to zoom 14 (sometimes 15), and then the browser/app client uses the vector data to scale up at subsequent zoom levels.


=====================================
include/coordinates.h
=====================================
@@ -9,8 +9,10 @@
 
 #ifdef FAT_TILE_INDEX
 typedef uint32_t TileCoordinate;
+#define TILE_COORDINATE_MAX UINT32_MAX
 #else
 typedef uint16_t TileCoordinate;
+#define TILE_COORDINATE_MAX UINT16_MAX
 #endif
 class TileCoordinates_ {
 


=====================================
include/geom.h
=====================================
@@ -89,5 +89,10 @@ void make_valid(GeometryT &geom) { }
 
 void make_valid(MultiPolygon &mp);
 
+Point intersect_edge(Point const &a, Point const &b, char edge, Box const &bbox);
+char bit_code(Point const &p, Box const &bbox);
+void fast_clip(Ring &points, Box const &bbox);
+void fast_clip(MultiPolygon &mp, Box const &bbox);
+
 #endif //_GEOM_TYPES_H
 


=====================================
include/osm_lua_processing.h
=====================================
@@ -161,6 +161,7 @@ public:
 	void AttributeBooleanWithMinZoom(const std::string &key, const bool val, const char minzoom);
 	void MinZoom(const double z);
 	void ZOrder(const double z);
+	void ZOrderWithScale(const double z, const double scale);
 	
 	// Relation scan support
 	kaguya::optional<int> NextRelation();


=====================================
include/output_object.h
=====================================
@@ -15,6 +15,12 @@
 #include "osmformat.pb.h"
 #include "vector_tile.pb.h"
 
+#ifdef FLOAT_Z_ORDER
+typedef float ZOrder;
+#else
+typedef int8_t ZOrder;
+#endif
+
 enum OutputGeometryType : unsigned int { POINT_, LINESTRING_, MULTILINESTRING_, POLYGON_ };
 
 #define OSMID_TYPE_OFFSET	40
@@ -46,16 +52,18 @@ protected:
 public:
 	NodeID objectID 			: 42;					// id of way (linestring/polygon) or node (point)
 	uint_least8_t layer 		: 8;					// what layer is it in?
-	int8_t z_order				: 8;					// z_order: used for sorting features within layers
+	ZOrder z_order				;						// z_order: used for sorting features within layers
 	OutputGeometryType geomType : 2;					// point, linestring, polygon
 	unsigned minZoom 			: 4;
 
 	AttributeStoreRef attributes;
 
-	void setZOrder(const int z) {
+	void setZOrder(const ZOrder z) {
+#ifndef FLOAT_Z_ORDER
 		if (z <= -127 || z >= 127) {
 			throw std::runtime_error("z_order is limited to 1 byte signed integer.");
 		}
+#endif
 		z_order = z;
 	}
 
@@ -156,14 +164,6 @@ LatpLon buildNodeGeometry(OSMStore &osmStore, OutputObject const &oo, const Tile
 
 bool operator==(const OutputObjectRef x, const OutputObjectRef y);
 
-/**
- * Do lexicographic comparison, with the order of: layer, geomType, attributes, and objectID.
- * Note that attributes is preferred to objectID.
- * It is to arrange objects with the identical attributes continuously.
- * Such objects will be merged into one object, to reduce the size of output.
- */
-bool operator<(const OutputObjectRef x, const OutputObjectRef y);
-
 namespace vector_tile {
 	bool operator==(const vector_tile::Tile_Value &x, const vector_tile::Tile_Value &y);
 	bool operator<(const vector_tile::Tile_Value &x, const vector_tile::Tile_Value &y);


=====================================
include/shared_data.h
=====================================
@@ -24,6 +24,7 @@ struct LayerDef {
 	uint filterBelow;
 	double filterArea;
 	uint combinePolygonsBelow;
+	bool sortZOrderAscending;
 	std::string source;
 	std::vector<std::string> sourceColumns;
 	bool allSourceColumns;
@@ -44,14 +45,14 @@ public:
 	// Define a layer (as read from the .json file)
 	uint addLayer(std::string name, uint minzoom, uint maxzoom,
 			uint simplifyBelow, double simplifyLevel, double simplifyLength, double simplifyRatio, 
-			uint filterBelow, double filterArea, uint combinePolygonsBelow,
+			uint filterBelow, double filterArea, uint combinePolygonsBelow, bool sortZOrderAscending,
 			const std::string &source,
 			const std::vector<std::string> &sourceColumns,
 			bool allSourceColumns,
 			bool indexed,
 			const std::string &indexName,
 			const std::string &writeTo);
-
+	std::vector<bool> getSortOrders();
 	rapidjson::Value serialiseToJSONValue(rapidjson::Document::AllocatorType &allocator) const;
 	std::string serialiseToJSON() const;
 };


=====================================
include/shp_mem_tiles.h
=====================================
@@ -14,7 +14,7 @@ public:
 	void CreateNamedLayerIndex(const std::string &layerName);
 
 	// Used in shape file loading
-	OutputObjectRef AddObject(uint_least8_t layerNum,
+	OutputObjectRef StoreShapefileGeometry(uint_least8_t layerNum,
 		const std::string &layerName, 
 		enum OutputGeometryType geomType,
 		Geometry geometry, 
@@ -23,6 +23,10 @@ public:
 	void AddObject(TileCoordinates const &index, OutputObjectRef const &oo) {
 		tileIndex[index].push_back(oo);
 	}
+	void AddObjectToLargeIndex(Box const &envelope, OutputObjectRef const &oo) {
+		std::lock_guard<std::mutex> lock(mutex);
+		box_rtree.insert(std::make_pair(envelope, oo));
+	}
 	std::vector<uint> QueryMatchingGeometries(const std::string &layerName, bool once, Box &box, 
 		std::function<std::vector<IndexValue>(const RTree &rtree)> indexQuery, 
 		std::function<bool(OutputObject const &oo)> checkQuery) const;


=====================================
include/tile_data.h
=====================================
@@ -19,6 +19,10 @@ protected:
 	std::mutex mutex;
 	TileIndex tileIndex;
 	std::deque<OutputObject> objects;
+	
+	// rtree index of large objects
+	using oo_rtree_param_type = boost::geometry::index::quadratic<128>;
+	boost::geometry::index::rtree< std::pair<Box,OutputObjectRef>, oo_rtree_param_type> box_rtree;
 
 	unsigned int baseZoom;
 
@@ -32,6 +36,8 @@ public:
 		MergeTileCoordsAtZoom(zoom, baseZoom, tileIndex, dstCoords);
 	}
 
+	void MergeLargeCoordsAtZoom(uint zoom, TileCoordinatesSet &dstCoords);
+
 	///This must be thread safe!
 	void MergeSingleTileDataAtZoom(TileCoordinates dstIndex, uint zoom, std::vector<OutputObjectRef> &dstTile) {
 		MergeSingleTileDataAtZoom(dstIndex, zoom, baseZoom, tileIndex, dstTile);
@@ -48,6 +54,13 @@ public:
 		tileIndex[index].push_back(oo);
 	}
 
+	void AddObjectToLargeIndex(Box const &envelope, OutputObjectRef const &oo) {
+		std::lock_guard<std::mutex> lock(mutex);
+		box_rtree.insert(std::make_pair(envelope, oo));
+	}
+
+	void MergeLargeObjects(TileCoordinates dstIndex, uint zoom, std::vector<OutputObjectRef> &dstTile);
+
 private:	
 	static void MergeTileCoordsAtZoom(uint zoom, uint baseZoom, const TileIndex &srcTiles, TileCoordinatesSet &dstCoords);
 	static void MergeSingleTileDataAtZoom(TileCoordinates dstIndex, uint zoom, uint baseZoom, const TileIndex &srcTiles, std::vector<OutputObjectRef> &dstTile);
@@ -55,7 +68,9 @@ private:
 
 TileCoordinatesSet GetTileCoordinates(std::vector<class TileDataSource *> const &sources, unsigned int zoom);
 
-std::vector<OutputObjectRef> GetTileData(std::vector<class TileDataSource *> const &sources, TileCoordinates coordinates, unsigned int zoom);
+std::vector<OutputObjectRef> GetTileData(std::vector<class TileDataSource *> const &sources,
+                                         std::vector<bool> const &sortOrders, 
+                                         TileCoordinates coordinates, unsigned int zoom);
 
 OutputObjectsConstItPair GetObjectsAtSubLayer(std::vector<OutputObjectRef> const &data, uint_least8_t layerNum);
 


=====================================
resources/process-openmaptiles.lua
=====================================
@@ -133,10 +133,12 @@ majorRoadValues = Set { "motorway", "trunk", "primary" }
 mainRoadValues  = Set { "secondary", "motorway_link", "trunk_link", "primary_link", "secondary_link" }
 midRoadValues   = Set { "tertiary", "tertiary_link" }
 minorRoadValues = Set { "unclassified", "residential", "road", "living_street" }
-trackValues     = Set { "cycleway", "byway", "bridleway", "track" }
-pathValues      = Set { "footway", "path", "steps", "pedestrian" }
+trackValues     = Set { "track" }
+pathValues      = Set { "footway", "cycleway", "bridleway", "path", "steps", "pedestrian" }
 linkValues      = Set { "motorway_link", "trunk_link", "primary_link", "secondary_link", "tertiary_link" }
 constructionValues = Set { "primary", "secondary", "tertiary", "motorway", "service", "trunk", "track" }
+pavedValues     = Set { "paved", "asphalt", "cobblestone", "concrete", "concrete:lanes", "concrete:plates", "metal", "paving_stones", "sett", "unhewn_cobblestone", "wood" }
+unpavedValues   = Set { "unpaved", "compacted", "dirt", "earth", "fine_gravel", "grass", "grass_paver", "gravel", "gravel_turf", "ground", "ice", "mud", "pebblestone", "salt", "sand", "snow", "woodchips" }
 
 aerowayBuildings= Set { "terminal", "gate", "tower" }
 landuseKeys     = Set { "school", "university", "kindergarten", "college", "library", "hospital",
@@ -285,7 +287,7 @@ function way_function(way)
 	-- Roads ('transportation' and 'transportation_name', plus 'transportation_name_detail')
 	if highway~="" then
 		local access = way:Find("access")
-		if access=="private" or access=="no" then return end
+		local surface = way:Find("surface")
 
 		local h = highway
 		local minzoom = 99
@@ -332,6 +334,9 @@ function way_function(way)
 			way:Attribute("class", h)
 			SetBrunnelAttributes(way)
 			if ramp then way:AttributeNumeric("ramp",1) end
+			if access=="private" or access=="no" then way:Attribute("access", "no") end
+			if pavedValues[surface] then way:Attribute("surface", "paved") end
+			if unpavedValues[surface] then way:Attribute("surface", "unpaved") end
 
 			-- Service
 			if highway == "service" and service ~="" then way:Attribute("service", service) end


=====================================
src/geom.cpp
=====================================
@@ -142,3 +142,80 @@ void make_valid(MultiPolygon &mp)
 	}
 	mp = result;
 }
+
+// ---------------
+// Sutherland-Hodgeman clipper
+// ported from Volodymyr Agafonkin's https://github.com/mapbox/lineclip
+
+// intersect a segment against one of the 4 lines that make up the bbox
+Point intersect_edge(Point const &a, Point const &b, char edge, Box const &bbox) {
+	if (edge & 8) return Point(a.x() + (b.x() - a.x()) * (bbox.max_corner().y() - a.y()) / (b.y() - a.y()), bbox.max_corner().y()); // top
+	if (edge & 4) return Point(a.x() + (b.x() - a.x()) * (bbox.min_corner().y() - a.y()) / (b.y() - a.y()), bbox.min_corner().y()); // bottom
+	if (edge & 2) return Point(bbox.max_corner().x(), a.y() + (b.y() - a.y()) * (bbox.max_corner().x() - a.x()) / (b.x() - a.x())); // right
+	if (edge & 1) return Point(bbox.min_corner().x(), a.y() + (b.y() - a.y()) * (bbox.min_corner().x() - a.x()) / (b.x() - a.x())); // left
+	throw std::runtime_error("intersect called on non-intersection");
+}
+
+// bit code reflects the point position relative to the bbox:
+//         left  mid  right
+//    top  1001  1000  1010
+//    mid  0001  0000  0010
+// bottom  0101  0100  0110
+
+char bit_code(Point const &p, Box const &bbox) {
+	char code = 0;
+	if      (p.x() < bbox.min_corner().x()) code |= 1; // left
+	else if (p.x() > bbox.max_corner().x()) code |= 2; // right
+	if      (p.y() < bbox.min_corner().y()) code |= 4; // bottom
+	else if (p.y() > bbox.max_corner().y()) code |= 8; // top
+	return code;
+}
+
+// Sutherland-Hodgeman polygon clipping algorithm
+void fast_clip(Ring &points, Box const &bbox) {
+	// clip against each side of the clip rectangle
+	for (char edge = 1; edge <= 8; edge *= 2) {
+		Ring result;
+		Point prev = points[points.size() - 1];
+		bool prevInside = (bit_code(prev, bbox) & edge)==0;
+
+		for (unsigned int i = 0; i<points.size(); i++) {
+			Point p = points[i];
+			bool inside = (bit_code(p, bbox) & edge)==0;
+
+			// if segment goes through the clip window, add an intersection
+			if (inside!=prevInside) result.emplace_back(intersect_edge(prev, p, edge, bbox));
+			if (inside) result.emplace_back(p); // add a point if it's inside
+			prev = p;
+			prevInside = inside;
+		}
+		points = std::move(result);
+		if (points.size()==0) break;
+	}
+}
+
+// Wrappers for polygon/multipolygon
+void fast_clip(Polygon &polygon, Box const &bbox) {
+	fast_clip(polygon.outer(), bbox);
+	if (polygon.outer().empty()) {
+		polygon.inners().resize(0);
+		return;
+	}
+	for (auto &inner: polygon.inners()) {
+		fast_clip(inner, bbox);
+	}
+	polygon.inners().erase(std::remove_if(
+		polygon.inners().begin(), polygon.inners().end(), 
+		[](const Ring &inner) -> bool { return inner.empty(); }),
+		polygon.inners().end());
+}
+
+void fast_clip(MultiPolygon &mp, Box const &bbox) {
+	for (auto &polygon: mp) {
+		fast_clip(polygon, bbox);
+	}
+	mp.erase(std::remove_if(
+		mp.begin(), mp.end(), 
+		[](const Polygon &poly) -> bool { return poly.outer().empty(); }),
+		mp.end());
+}


=====================================
src/osm_lua_processing.cpp
=====================================
@@ -55,7 +55,7 @@ OsmLuaProcessing::OsmLuaProcessing(
 		.addOverloadedFunctions("AttributeNumeric", &OsmLuaProcessing::AttributeNumeric, &OsmLuaProcessing::AttributeNumericWithMinZoom)
 		.addOverloadedFunctions("AttributeBoolean", &OsmLuaProcessing::AttributeBoolean, &OsmLuaProcessing::AttributeBooleanWithMinZoom)
 		.addFunction("MinZoom", &OsmLuaProcessing::MinZoom)
-		.addFunction("ZOrder", &OsmLuaProcessing::ZOrder)
+		.addOverloadedFunctions("ZOrder", &OsmLuaProcessing::ZOrder, &OsmLuaProcessing::ZOrderWithScale)
 		.addFunction("Accept", &OsmLuaProcessing::Accept)
 		.addFunction("NextRelation", &OsmLuaProcessing::NextRelation)
 		.addFunction("FindInRelation", &OsmLuaProcessing::FindInRelation)
@@ -501,7 +501,21 @@ void OsmLuaProcessing::MinZoom(const double z) {
 // Set z_order
 void OsmLuaProcessing::ZOrder(const double z) {
 	if (outputs.size()==0) { ProcessingError("Can't set z_order if no Layer set"); return; }
+#ifdef FLOAT_Z_ORDER
+	outputs.back().first->setZOrder(make_valid<float>(z));
+#else
 	outputs.back().first->setZOrder(make_valid<int>(z));
+#endif
+}
+
+// Set z_order (variant with scaling)
+void OsmLuaProcessing::ZOrderWithScale(const double z, const double scale) {
+	if (outputs.size()==0) { ProcessingError("Can't set z_order if no Layer set"); return; }
+#ifdef FLOAT_Z_ORDER
+	outputs.back().first->setZOrder(make_valid<float>(z));
+#else
+	outputs.back().first->setZOrder(make_valid<int>(z/scale*127));
+#endif
 }
 
 // Read scanned relations
@@ -609,29 +623,45 @@ void OsmLuaProcessing::setWay(WayID wayId, LatpLonVec const &llVec, const tag_ma
 		// create a list of tiles this way passes through (tileSet)
 		unordered_set<TileCoordinates> tileSet;
 		try {
-			insertIntermediateTiles(osmStore.llListLinestring(llVecPtr->cbegin(),llVecPtr->cend()), this->config.baseZoom, tileSet);
+			Linestring ls = osmStore.llListLinestring(llVecPtr->cbegin(),llVecPtr->cend());
+			insertIntermediateTiles(ls, this->config.baseZoom, tileSet);
 
 			// then, for each tile, store the OutputObject for each layer
 			bool polygonExists = false;
+			TileCoordinate minTileX = TILE_COORDINATE_MAX, maxTileX = 0, minTileY = TILE_COORDINATE_MAX, maxTileY = 0;
 			for (auto it = tileSet.begin(); it != tileSet.end(); ++it) {
 				TileCoordinates index = *it;
+				minTileX = std::min(index.x, minTileX);
+				minTileY = std::min(index.y, minTileY);
+				maxTileX = std::max(index.x, maxTileX);
+				maxTileY = std::max(index.y, maxTileY);
 				for (auto jt = this->outputs.begin(); jt != this->outputs.end(); ++jt) {
 					if (jt->first->geomType == POLYGON_) {
 						polygonExists = true;
 						continue;
 					}
-					osmMemTiles.AddObject(index, jt->first);
+					osmMemTiles.AddObject(index, jt->first); // not a polygon
 				}
 			}
 
 			// for polygon, fill inner tiles
 			if (polygonExists) {
-				fillCoveredTiles(tileSet);
-				for (auto it = tileSet.begin(); it != tileSet.end(); ++it) {
-					TileCoordinates index = *it;
-					for (auto jt = this->outputs.begin(); jt != this->outputs.end(); ++jt) {
-						if (jt->first->geomType != POLYGON_) continue;
-						osmMemTiles.AddObject(index, jt->first);
+				bool tilesetFilled = false;
+				uint size = (maxTileX - minTileX + 1) * (maxTileY - minTileY + 1);
+				for (auto jt = this->outputs.begin(); jt != this->outputs.end(); ++jt) {
+					if (jt->first->geomType != POLYGON_) continue;
+					if (size>= 16) {
+						// Larger objects - add to rtree
+						Box box = Box(geom::make<Point>(minTileX, minTileY),
+						              geom::make<Point>(maxTileX, maxTileY));
+						osmMemTiles.AddObjectToLargeIndex(box, jt->first);
+					} else {
+						// Smaller objects - add to each individual tile index
+						if (!tilesetFilled) { fillCoveredTiles(tileSet); tilesetFilled = true; }
+						for (auto it = tileSet.begin(); it != tileSet.end(); ++it) {
+							TileCoordinates index = *it;
+							osmMemTiles.AddObject(index, jt->first);
+						}
 					}
 				}
 			}
@@ -676,28 +706,46 @@ void OsmLuaProcessing::setRelation(int64_t relationId, WayVec const &outerWayVec
 			return;
 		}		
 
-		unordered_set<TileCoordinates> tileSet;
-		if (mp.size() == 1) {
-			insertIntermediateTiles(mp[0].outer(), this->config.baseZoom, tileSet);
-			fillCoveredTiles(tileSet);
-		} else {
-			for (Polygon poly: mp) {
-				unordered_set<TileCoordinates> tileSetTmp;
-				insertIntermediateTiles(poly.outer(), this->config.baseZoom, tileSetTmp);
-				fillCoveredTiles(tileSetTmp);
-				tileSet.insert(tileSetTmp.begin(), tileSetTmp.end());
-			}
-		}
-
 		for (auto jt = this->outputs.begin(); jt != this->outputs.end(); ++jt) {
 			// Store the attributes of the generated geometry
 			jt->first->setAttributeSet(attributeStore.store_set(jt->second));		
 		}
 
+		unordered_set<TileCoordinates> tileSet;
+		bool singleOuter = mp.size()==1;
+		for (Polygon poly: mp) {
+			unordered_set<TileCoordinates> tileSetTmp;
+			insertIntermediateTiles(poly.outer(), this->config.baseZoom, tileSetTmp);
+			fillCoveredTiles(tileSetTmp);
+			if (singleOuter) {
+				tileSet = std::move(tileSetTmp);
+			} else {
+				tileSet.insert(tileSetTmp.begin(), tileSetTmp.end());
+			}
+		}
+		
+		TileCoordinate minTileX = TILE_COORDINATE_MAX, maxTileX = 0, minTileY = TILE_COORDINATE_MAX, maxTileY = 0;
 		for (auto it = tileSet.begin(); it != tileSet.end(); ++it) {
 			TileCoordinates index = *it;
-			for (auto jt = this->outputs.begin(); jt != this->outputs.end(); ++jt) {
-				osmMemTiles.AddObject(index, jt->first);
+			minTileX = std::min(index.x, minTileX);
+			minTileY = std::min(index.y, minTileY);
+			maxTileX = std::max(index.x, maxTileX);
+			maxTileY = std::max(index.y, maxTileY);
+		}
+		for (auto jt = this->outputs.begin(); jt != this->outputs.end(); ++jt) {
+			if (tileSet.size()>=16) {
+				// Larger objects - add to rtree
+				// note that the bbox is currently the envelope of the entire multipolygon,
+				// which is suboptimal in shapes like (_) ...... (_) where the outers are significantly disjoint
+				Box box = Box(geom::make<Point>(minTileX, minTileY),
+				              geom::make<Point>(maxTileX, maxTileY));
+				osmMemTiles.AddObjectToLargeIndex(box, jt->first);
+			} else {
+				// Smaller objects - add to each individual tile index
+				for (auto it = tileSet.begin(); it != tileSet.end(); ++it) {
+					TileCoordinates index = *it;
+					osmMemTiles.AddObject(index, jt->first);
+				}
 			}
 		}
 


=====================================
src/output_object.cpp
=====================================
@@ -167,23 +167,11 @@ Geometry buildWayGeometry(OSMStore &osmStore, OutputObject const &oo, const Tile
 					std::min(box.max_corner().y(), extBox.max_corner().y()));
 			}
 
-			Polygon clippingPolygon;
-			geom::convert(box, clippingPolygon);
-	
-			try {
-				MultiPolygon mp;
-				if (geom::within(input, clippingPolygon)) { geom::assign(mp, input); return mp; }
-				// Work around boost::geometry intersection issues by doing each constituent polygon at a time
-				for (auto &p : input) {
-					MultiPolygon out;
-					geom::intersection(p, clippingPolygon, out);
-					for (auto &o : out) mp.emplace_back(move(o));
-				}
-				return mp;
-			} catch (geom::overlay_invalid_input_exception &err) {
-				std::cout << "Couldn't clip polygon (self-intersection)" << std::endl;
-				return MultiPolygon(); // blank
-			}
+			MultiPolygon mp;
+			geom::assign(mp, input);
+			fast_clip(mp, box);
+			geom::correct(mp);
+			return mp;
 		}
 
 		default:
@@ -234,23 +222,6 @@ bool operator==(const OutputObjectRef x, const OutputObjectRef y) {
 		x->objectID == y->objectID;
 } 
 
-// Do lexicographic comparison, with the order of: layer, geomType, attributes, and objectID.
-// Note that attributes is preffered to objectID.
-// It is to arrange objects with the identical attributes continuously.
-// Such objects will be merged into one object, to reduce the size of output.
-bool operator<(const OutputObjectRef x, const OutputObjectRef y) {
-	if (x->layer < y->layer) return true;
-	if (x->layer > y->layer) return false;
-	if (x->z_order < y->z_order) return true;
-	if (x->z_order > y->z_order) return false;
-	if (x->geomType < y->geomType) return true;
-	if (x->geomType > y->geomType) return false;
-	if (x->attributes.get() < y->attributes.get()) return true;
-	if (x->attributes.get() > y->attributes.get()) return false;
-	if (x->objectID < y->objectID) return true;
-	return false;
-} 
-
 namespace vector_tile {
 	bool operator==(const vector_tile::Tile_Value &x, const vector_tile::Tile_Value &y) {
 		std::string strx = x.SerializeAsString();


=====================================
src/read_shp.cpp
=====================================
@@ -177,7 +177,7 @@ void readShapefile(const Box &clippingBox,
 				if (indexField>-1) { name=DBFReadStringAttribute(dbf, i, indexField); hasName = true;}
 
 				auto &attributeStore = osmLuaProcessing.getAttributeStore();
-				OutputObjectRef oo = shpMemTiles.AddObject(layerNum, layerName, POINT_, p, isIndexed, hasName, name, attributeStore.empty_set(), layer.minzoom);
+				OutputObjectRef oo = shpMemTiles.StoreShapefileGeometry(layerNum, layerName, POINT_, p, isIndexed, hasName, name, attributeStore.empty_set(), layer.minzoom);
 
 				addShapefileAttributes(dbf, oo, i, columnMap, columnTypeMap, layers, osmLuaProcessing);
 			}
@@ -199,7 +199,7 @@ void readShapefile(const Box &clippingBox,
 					if (indexField>-1) { name=DBFReadStringAttribute(dbf, i, indexField); hasName = true;}
 
 					auto &attributeStore = osmLuaProcessing.getAttributeStore();
-					OutputObjectRef oo = shpMemTiles.AddObject(layerNum, layerName, LINESTRING_, *it, isIndexed, hasName, name, attributeStore.empty_set(), layer.minzoom);
+					OutputObjectRef oo = shpMemTiles.StoreShapefileGeometry(layerNum, layerName, LINESTRING_, *it, isIndexed, hasName, name, attributeStore.empty_set(), layer.minzoom);
 
 					addShapefileAttributes(dbf, oo, i, columnMap, columnTypeMap, layers, osmLuaProcessing);
 				}
@@ -272,7 +272,7 @@ void readShapefile(const Box &clippingBox,
 
 				// create OutputObject
 				auto &attributeStore = osmLuaProcessing.getAttributeStore();
-				OutputObjectRef oo = shpMemTiles.AddObject(layerNum, layerName, POLYGON_, out, isIndexed, hasName, name, attributeStore.empty_set(), layer.minzoom);
+				OutputObjectRef oo = shpMemTiles.StoreShapefileGeometry(layerNum, layerName, POLYGON_, out, isIndexed, hasName, name, attributeStore.empty_set(), layer.minzoom);
 
 				addShapefileAttributes(dbf, oo, i, columnMap, columnTypeMap, layers, osmLuaProcessing);
 			}


=====================================
src/shared_data.cpp
=====================================
@@ -19,7 +19,7 @@ SharedData::~SharedData() { }
 // Define a layer (as read from the .json file)
 uint LayerDefinition::addLayer(string name, uint minzoom, uint maxzoom,
 		uint simplifyBelow, double simplifyLevel, double simplifyLength, double simplifyRatio, 
-		uint filterBelow, double filterArea, uint combinePolygonsBelow,
+		uint filterBelow, double filterArea, uint combinePolygonsBelow, bool sortZOrderAscending,
 		const std::string &source,
 		const std::vector<std::string> &sourceColumns,
 		bool allSourceColumns,
@@ -29,7 +29,7 @@ uint LayerDefinition::addLayer(string name, uint minzoom, uint maxzoom,
 
 	bool isWriteTo = !writeTo.empty();
 	LayerDef layer = { name, minzoom, maxzoom, simplifyBelow, simplifyLevel, simplifyLength, simplifyRatio, 
-		filterBelow, filterArea, combinePolygonsBelow,
+		filterBelow, filterArea, combinePolygonsBelow, sortZOrderAscending,
 		source, sourceColumns, allSourceColumns, indexed, indexName,
 		std::map<std::string,uint>(), isWriteTo };
 	layers.push_back(layer);
@@ -54,6 +54,11 @@ uint LayerDefinition::addLayer(string name, uint minzoom, uint maxzoom,
 	return layerNum;
 }
 
+std::vector<bool> LayerDefinition::getSortOrders() {
+	std::vector<bool> orders;
+	for (auto &layer : layers) { orders.emplace_back(layer.sortZOrderAscending); }
+	return orders;
+}
 
 Value LayerDefinition::serialiseToJSONValue(rapidjson::Document::AllocatorType &allocator) const {
 	Value layerArray(kArrayType);
@@ -129,6 +134,13 @@ void Config::readConfig(rapidjson::Document &jsonConfig, bool &hasClippingBox, B
 		cerr << "\"compress\" should be any of \"gzip\",\"deflate\",\"none\" in JSON file." << endl;
 		exit (EXIT_FAILURE);
 	}
+	if (endZoom>15) {
+		cout << "**** WARNING ****" << endl;
+		cout << "You're generating tiles up to z" << endZoom << ". You probably don't want to do that." << endl;
+		cout << "Standard practice is to generate vector tiles up to z14. Your renderer will 'overzoom' the z14 tiles for higher resolutions." << endl;
+		cout << "tilemaker may have excessive memory, time, and space requirements at higher zooms. You can find more information in the docs/ folder." << endl;
+		cout << "**** WARNING ****" << endl;
+	}
 #ifndef FAT_TILE_INDEX
 	if (endZoom>16) {
 		cerr << "Compile tilemaker with -DFAT_TILE_INDEX to enable tile output at zoom level 17 or greater" << endl;
@@ -191,6 +203,7 @@ void Config::readConfig(rapidjson::Document &jsonConfig, bool &hasClippingBox, B
 		int    filterBelow    = it->value.HasMember("filter_below"   ) ? it->value["filter_below"   ].GetInt()    : 0;
 		double filterArea     = it->value.HasMember("filter_area"    ) ? it->value["filter_area"    ].GetDouble() : 0.5;
 		int    combinePolyBelow=it->value.HasMember("combine_polygons_below") ? it->value["combine_polygons_below"].GetInt() : 0;
+		bool sortZOrderAscending = it->value.HasMember("z_order_ascending") ? it->value["z_order_ascending"].GetBool() : true;
 		string source = it->value.HasMember("source") ? it->value["source"].GetString() : "";
 		vector<string> sourceColumns;
 		bool allSourceColumns = false;
@@ -210,7 +223,7 @@ void Config::readConfig(rapidjson::Document &jsonConfig, bool &hasClippingBox, B
 
 		layers.addLayer(layerName, minZoom, maxZoom,
 				simplifyBelow, simplifyLevel, simplifyLength, simplifyRatio, 
-				filterBelow, filterArea, combinePolyBelow,
+				filterBelow, filterArea, combinePolyBelow, sortZOrderAscending,
 				source, sourceColumns, allSourceColumns, indexed, indexName,
 				writeTo);
 


=====================================
src/shp_mem_tiles.cpp
=====================================
@@ -51,7 +51,7 @@ void ShpMemTiles::CreateNamedLayerIndex(const std::string &layerName) {
 	indices[layerName]=RTree();
 }
 
-OutputObjectRef ShpMemTiles::AddObject(uint_least8_t layerNum,
+OutputObjectRef ShpMemTiles::StoreShapefileGeometry(uint_least8_t layerNum,
 	const std::string &layerName, enum OutputGeometryType geomType,
 	Geometry geometry, bool isIndexed, bool hasName, const std::string &name, AttributeStoreRef attributes, uint minzoom) {
 
@@ -116,15 +116,23 @@ OutputObjectRef ShpMemTiles::AddObject(uint_least8_t layerNum,
 }
 
 // Add an OutputObject to all tiles between min/max lat/lon
+// (only used for polygons)
 void ShpMemTiles::addToTileIndexByBbox(OutputObjectRef &oo, double minLon, double minLatp, double maxLon, double maxLatp) {
 	uint minTileX =  lon2tilex(minLon, baseZoom);
 	uint maxTileX =  lon2tilex(maxLon, baseZoom);
 	uint minTileY = latp2tiley(minLatp, baseZoom);
 	uint maxTileY = latp2tiley(maxLatp, baseZoom);
-	for (uint x=min(minTileX,maxTileX); x<=max(minTileX,maxTileX); x++) {
-		for (uint y=min(minTileY,maxTileY); y<=max(minTileY,maxTileY); y++) {
-			TileCoordinates index(x, y);
-			AddObject(index, oo);
+	uint size = (maxTileX - minTileX + 1) * (minTileY - maxTileY + 1);
+	if (size>=16) { 
+		// Larger objects - add to rtree
+		AddObjectToLargeIndex(Box(Point(minTileX, maxTileY), Point(maxTileX, minTileY)), oo);
+	} else {
+		// Smaller objects - add to each individual tile index
+		for (uint x=min(minTileX,maxTileX); x<=max(minTileX,maxTileX); x++) {
+			for (uint y=min(minTileY,maxTileY); y<=max(minTileY,maxTileY); y++) {
+				TileCoordinates index(x, y);
+				AddObject(index, oo);
+			}
 		}
 	}
 }


=====================================
src/tile_data.cpp
=====================================
@@ -29,6 +29,24 @@ void TileDataSource::MergeTileCoordsAtZoom(uint zoom, uint baseZoom, const TileI
 	}
 }
 
+// Find the tiles used by the "large objects" from the rtree index
+void TileDataSource::MergeLargeCoordsAtZoom(uint zoom, TileCoordinatesSet &dstCoords) {
+	for(auto const &result: box_rtree) {
+		int scale = pow(2, baseZoom-zoom);
+		TileCoordinate minx = result.first.min_corner().x() / scale;
+		TileCoordinate maxx = result.first.max_corner().x() / scale;
+		TileCoordinate miny = result.first.min_corner().y() / scale;
+		TileCoordinate maxy = result.first.max_corner().y() / scale;
+		for (int x=minx; x<=maxx; x++) {
+			for (int y=miny; y<=maxy; y++) {
+				TileCoordinates newIndex(x, y);
+				dstCoords.insert(newIndex);
+			}
+		}
+	}
+}
+
+// Copy objects from the tile at dstIndex (in the dataset srcTiles) into dstTile
 void TileDataSource::MergeSingleTileDataAtZoom(TileCoordinates dstIndex, uint zoom, uint baseZoom, const TileIndex &srcTiles, std::vector<OutputObjectRef> &dstTile) {
 	if (zoom==baseZoom) {
 		// at z14, we can just use tileIndex
@@ -59,24 +77,55 @@ void TileDataSource::MergeSingleTileDataAtZoom(TileCoordinates dstIndex, uint zo
 	}
 }
 
+// Copy objects from the large index into dstTile
+void TileDataSource::MergeLargeObjects(TileCoordinates dstIndex, uint zoom, std::vector<OutputObjectRef> &dstTile) {
+	int scale = pow(2, baseZoom - zoom);
+	TileCoordinates srcIndex1( dstIndex.x   *scale  ,  dstIndex.y   *scale  );
+	TileCoordinates srcIndex2((dstIndex.x+1)*scale-1, (dstIndex.y+1)*scale-1);
+	Box box = Box(geom::make<Point>(srcIndex1.x, srcIndex1.y),
+	              geom::make<Point>(srcIndex2.x, srcIndex2.y));
+	for(auto const &result: box_rtree | boost::geometry::index::adaptors::queried(boost::geometry::index::intersects(box)))
+		dstTile.push_back(result.second);
+}
+
 TileCoordinatesSet GetTileCoordinates(std::vector<class TileDataSource *> const &sources, unsigned int zoom) {
 	TileCoordinatesSet tileCoordinates;
 
 	// Create list of tiles
 	tileCoordinates.clear();
-	for(size_t i=0; i<sources.size(); i++)
+	for(size_t i=0; i<sources.size(); i++) {
 		sources[i]->MergeTileCoordsAtZoom(zoom, tileCoordinates);
+		sources[i]->MergeLargeCoordsAtZoom(zoom, tileCoordinates);
+	}
 
 	return tileCoordinates;
 }
 
-std::vector<OutputObjectRef> GetTileData(std::vector<class TileDataSource *> const &sources, TileCoordinates coordinates, unsigned int zoom)
-{
+std::vector<OutputObjectRef> GetTileData(std::vector<class TileDataSource *> const &sources, 
+                                         std::vector<bool> const &sortOrders, TileCoordinates coordinates, 
+                                         unsigned int zoom) {
 	std::vector<OutputObjectRef> data;
-	for(size_t i=0; i<sources.size(); i++)
+	for(size_t i=0; i<sources.size(); i++) {
 		sources[i]->MergeSingleTileDataAtZoom(coordinates, zoom, data);
+		sources[i]->MergeLargeObjects(coordinates, zoom, data);
+	}
 
-	boost::sort::pdqsort(data.begin(), data.end());
+	// Lexicographic comparison, with the order of: layer, geomType, attributes, and objectID.
+	// Note that attributes is preferred to objectID.
+	// It is to arrange objects with the identical attributes continuously.
+	// Such objects will be merged into one object, to reduce the size of output.
+	boost::sort::pdqsort(data.begin(), data.end(), [&sortOrders](const OutputObjectRef x, const OutputObjectRef y) -> bool {
+		if (x->layer < y->layer) return true;
+		if (x->layer > y->layer) return false;
+		if (x->z_order < y->z_order) return  sortOrders[x->layer];
+		if (x->z_order > y->z_order) return !sortOrders[x->layer];
+		if (x->geomType < y->geomType) return true;
+		if (x->geomType > y->geomType) return false;
+		if (x->attributes.get() < y->attributes.get()) return true;
+		if (x->attributes.get() > y->attributes.get()) return false;
+		if (x->objectID < y->objectID) return true;
+		return false;
+	});
 	data.erase(unique(data.begin(), data.end()), data.end());
 	return data;
 }


=====================================
src/tilemaker.cpp
=====================================
@@ -326,6 +326,7 @@ int main(int argc, char* argv[]) {
 	// ----	Read all PBFs
 	
 	PbfReader pbfReader(osmStore);
+	std::vector<bool> sortOrders = layers.getSortOrders();
 
 	if (!mapsplit) {
 		for (auto inputFile : inputFiles) {
@@ -469,7 +470,7 @@ int main(int argc, char* argv[]) {
 				for(std::size_t i = start_index; i < end_index; ++i) {
 					unsigned int zoom = tile_coordinates[i].first;
 					TileCoordinates coords = tile_coordinates[i].second;
-					outputProc(pool, sharedData, osmStore, GetTileData(sources, coords, zoom), coords, zoom);
+					outputProc(pool, sharedData, osmStore, GetTileData(sources, sortOrders, coords, zoom), coords, zoom);
 				}
 
 				const std::lock_guard<std::mutex> lock(io_mutex);



View it on GitLab: https://salsa.debian.org/debian-gis-team/tilemaker/-/commit/1ca89b0c54838f3d65a539f7f37bef9e3de5d020

-- 
View it on GitLab: https://salsa.debian.org/debian-gis-team/tilemaker/-/commit/1ca89b0c54838f3d65a539f7f37bef9e3de5d020
You're receiving this email because of your account on salsa.debian.org.


-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://alioth-lists.debian.net/pipermail/pkg-grass-devel/attachments/20230701/f3c079de/attachment-0001.htm>


More information about the Pkg-grass-devel mailing list