Speed up runtime by reducing allocation churn in tile generation

include/coordinates_geom.h

@@ -23,7 +23,9 @@ class TileBbox {
 	TileBbox(TileCoordinates i, uint z, bool h, bool e);
 
 	std::pair<int,int> scaleLatpLon(double latp, double lon) const;
-	std::vector<Point> scaleRing(Ring const &src) const;
+	void scaleRing(Ring &dst, Ring const &src) const;
+	Ring scaleRing(Ring const &src) const;
+	void scaleGeometry(MultiPolygon &dst, MultiPolygon const &src) const;
 	MultiPolygon scaleGeometry(MultiPolygon const &src) const;
 	std::pair<double, double> floorLatpLon(double latp, double lon) const;
 

include/geometry/correct.hpp

@@ -10,6 +10,7 @@
  * ----------------------------------------------------------------------------
  */
 
+#include <utility>
 #include <vector>
 #include <boost/geometry.hpp>
 #include <boost/geometry/geometries/point_xy.hpp>
@@ -25,7 +26,7 @@ namespace impl {
 template<typename C, typename T>
 static inline void result_combine(C &result, T &&new_element)
 {
-    result.push_back(new_element);
+    result.push_back(std::forward<T>(new_element));
 
    	for(std::size_t i = 0; i < result.size() - 1; ) {
         if(!boost::geometry::intersects(result[i], result.back())) {
@@ -96,6 +97,8 @@ static inline void dissolve_find_intersections(
 	if(ring.empty()) return;
 
 	boost::geometry::index::rtree<std::pair< boost::geometry::model::segment<point_t>, std::size_t >, boost::geometry::index::quadratic<16>> index;
+	std::vector<point_t> output;
+	output.reserve(2);
 
 	// Generate all by-pass intersections in the graph
 	// Generate a list of all by-pass intersections
@@ -112,7 +115,7 @@ static inline void dissolve_find_intersections(
 			auto const &line_2 = iter.first;
 			auto j = iter.second;
 
-			std::vector<point_t> output;
+			output.clear();
 			boost::geometry::intersection(line_1, line_2, output);
 
 			for(auto const &p: output) {
@@ -281,9 +284,11 @@ static inline std::vector<ring_t> correct(ring_t const &ring, boost::geometry::o
 	dissolve_find_intersections(new_ring, pseudo_vertices, start_keys);
 
 	if(start_keys.empty()) {
-		if(std::abs(boost::geometry::area(new_ring)) > remove_spike_min_area) 
-			return { new_ring };
-		else
+		if(std::abs(boost::geometry::area(new_ring)) > remove_spike_min_area) {
+			std::vector<ring_t> result;
+			result.push_back(std::move(new_ring));
+			return result;
+		} else
 			return { };
 	}
 

include/osm_lua_processing.h

@@ -276,6 +276,8 @@ class OsmLuaProcessing {
 	const inline Point getPoint() {
 		return Point(lon/10000000.0,latp/10000000.0);
 	}
+
+	double projectedPolygonArea(const Polygon &p);
 
 	OSMStore &osmStore;	// global OSM store
 
@@ -310,6 +312,7 @@ class OsmLuaProcessing {
 	bool multiLinestringInited;
 	MultiPolygon multiPolygonCache;
 	bool multiPolygonInited;
+	geom::model::polygon<DegPoint> areaPolygonCache;
 
 	NodeID lastStoredGeometryId;
 	OutputGeometryType lastStoredGeometryType;

include/osm_store.h

@@ -332,6 +332,7 @@ class OSMStore
 	template<class WayIt>
 	Polygon llListPolygon(WayIt begin, WayIt end) const {
 		Polygon poly;
+		poly.outer().reserve(end - begin);
 		fillPoints(poly.outer(), begin, end);
 		boost::geometry::correct(poly);
 		return poly;
@@ -341,6 +342,7 @@ class OSMStore
 	template<class WayIt>
 	Linestring llListLinestring(WayIt begin, WayIt end) const {
 		Linestring ls;
+		ls.reserve(end - begin);
 		fillPoints(ls, begin, end);
 		return ls;
 	}

include/sharded_way_store.h

@@ -14,6 +14,7 @@ class ShardedWayStore : public WayStore {
 	void reopen() override;
 	void batchStart() override;
 	std::vector<LatpLon> at(WayID wayid) const override;
+	void at(WayID wayid, std::vector<LatpLon>& output) const override;
 	bool requiresNodes() const override;
 	void insertLatpLons(std::vector<WayStore::ll_element_t> &newWays) override;
 	void insertNodes(const std::vector<std::pair<WayID, std::vector<NodeID>>>& newWays) override;

include/sorted_way_store.h

@@ -88,6 +88,7 @@ class SortedWayStore: public WayStore {
 	void reopen() override;
 	void batchStart() override;
 	std::vector<LatpLon> at(WayID wayid) const override;
+	void at(WayID wayid, std::vector<LatpLon>& output) const override;
 	bool requiresNodes() const override { return true; }
 	void insertLatpLons(std::vector<WayStore::ll_element_t> &newWays) override;
 	void insertNodes(const std::vector<std::pair<WayID, std::vector<NodeID>>>& newWays) override;
@@ -107,6 +108,7 @@ class SortedWayStore: public WayStore {
 	);
 
 	static std::vector<NodeID> decodeWay(uint16_t flags, const uint8_t* input);
+	static void decodeWay(uint16_t flags, const uint8_t* input, std::vector<NodeID>& output);
 
 private:
 	bool compressWays;

include/tile_data.h

@@ -114,7 +114,7 @@ inline OutputObjectID outputObjectWithId<OutputObjectXYID>(const OutputObjectXYI
 
 template<typename OO> void collectLowZoomObjectsForTile(
 	const unsigned int& indexZoom,
-	typename std::vector<std::vector<OO>> objects,
+	const typename std::vector<std::vector<OO>>& objects,
 	unsigned int zoom,
 	const TileCoordinates& dstIndex,
 	std::vector<OutputObjectID>& output

include/way_store.h

@@ -15,6 +15,7 @@ class WayStore {
 	// meaningful for SortedWayStore
 	virtual void batchStart() = 0;
 	virtual std::vector<LatpLon> at(WayID wayid) const = 0;
+	virtual void at(WayID wayid, std::vector<LatpLon>& output) const = 0;
 	virtual bool requiresNodes() const = 0;
 	virtual void insertLatpLons(std::vector<ll_element_t>& newWays) = 0;
 	virtual void insertNodes(const std::vector<std::pair<WayID, std::vector<NodeID>>>& newWays) = 0;

include/way_stores.h

@@ -15,6 +15,7 @@ class BinarySearchWayStore: public WayStore {
 	void reopen() override;
 	void batchStart() override {}
 	std::vector<LatpLon> at(WayID wayid) const override;
+	void at(WayID wayid, std::vector<LatpLon>& output) const override;
 	bool requiresNodes() const override { return false; }
 	void insertLatpLons(std::vector<WayStore::ll_element_t> &newWays) override;
 	void insertNodes(const std::vector<std::pair<WayID, std::vector<NodeID>>>& newWays) override;

src/coordinates_geom.cpp

@@ -33,10 +33,10 @@ pair<int,int> TileBbox::scaleLatpLon(double latp, double lon) const {
 
 // Scaling with naive self-intersection check - if we've added the new point
 // within the last 5 points, then backtrack to the last time we added it
-std::vector<Point> TileBbox::scaleRing(Ring const &src) const {
-	std::vector<Point> points;
+void TileBbox::scaleRing(Ring &points, Ring const &src) const {
+	points.clear();
 	points.reserve(src.size());
-	for(auto &i: src) {
+	for(auto const &i: src) {
 		auto scaled = scaleLatpLon(i.y(), i.x()); // -> .first is x, .second is y
 		bool found = false;
 		for (size_t j=1; j<5; j++) {
@@ -48,36 +48,46 @@ std::vector<Point> TileBbox::scaleRing(Ring const &src) const {
 		}
 		if (!found) points.push_back(Point(scaled.first,scaled.second));
 	}
+}
+
+Ring TileBbox::scaleRing(Ring const &src) const {
+	Ring points;
+	scaleRing(points, src);
 	return points;
 }
 
-MultiPolygon TileBbox::scaleGeometry(MultiPolygon const &src) const {
-	MultiPolygon dst;
-	for(auto poly: src) {
-		Polygon p;
+void TileBbox::scaleGeometry(MultiPolygon &dst, MultiPolygon const &src) const {
+	if (dst.size() < src.size())
+		dst.resize(src.size());
+
+	size_t polygonCount = 0;
+	for(auto const &poly: src) {
+		Polygon &p = dst[polygonCount];
 
 		// Copy the outer ring
-		std::vector<Point> points = scaleRing(poly.outer());
-		if (points.size()<4) continue;
-		Ring outer;
-		geom::append(outer,points);
-		geom::append(p,outer);
+		scaleRing(p.outer(), poly.outer());
+		if (p.outer().size()<4)
+			continue;
 
 		// Copy the inner rings
-		int num_rings = 0;
-		for(auto &r: poly.inners()) {
-			points = scaleRing(r);
-			if (points.size()<4) continue;
-			Ring inner;
-			geom::append(inner,points);
-			num_rings++;
-			geom::interior_rings(p).resize(num_rings);
-			geom::append(p, inner, num_rings-1);
+		if (p.inners().size() < poly.inners().size())
+			p.inners().resize(poly.inners().size());
+		size_t innerCount = 0;
+		for(auto const &r: poly.inners()) {
+			Ring &points = p.inners()[innerCount];
+			scaleRing(points, r);
+			if (points.size()>=4)
+				innerCount++;
 		}
-
-		// Add to multipolygon
-		dst.push_back(p);
+		p.inners().resize(innerCount);
+		polygonCount++;
 	}
+	dst.resize(polygonCount);
+}
+
+MultiPolygon TileBbox::scaleGeometry(MultiPolygon const &src) const {
+	MultiPolygon dst;
+	scaleGeometry(dst, src);
 	return dst;
 }
 

src/geom.cpp

@@ -197,10 +197,20 @@ char bit_code(Point const &p, Box const &bbox) {
 }
 
 // Sutherland-Hodgeman polygon clipping algorithm
-void fast_clip(Ring &points, Box const &bbox) {
+void fast_clip(Ring &points, Box const &bbox, Ring &result) {
 	// clip against each side of the clip rectangle
 	for (char edge = 1; edge <= 8; edge *= 2) {
-		Ring result;
+		bool needsClip = false;
+		for (auto const &p: points) {
+			if (bit_code(p, bbox) & edge) {
+				needsClip = true;
+				break;
+			}
+		}
+		if (!needsClip) continue;
+
+		result.clear();
+		result.reserve(points.size() + 4);
 		Point prev = points[points.size() - 1];
 		bool prevInside = (bit_code(prev, bbox) & edge)==0;
 
@@ -214,20 +224,26 @@ void fast_clip(Ring &points, Box const &bbox) {
 			prev = p;
 			prevInside = inside;
 		}
-		points = std::move(result);
+		points.swap(result);
 		if (points.size()==0) break;
 	}
 }
 
+void fast_clip(Ring &points, Box const &bbox) {
+	Ring result;
+	fast_clip(points, bbox, result);
+}
+
 // Wrappers for polygon/multipolygon
 void fast_clip(Polygon &polygon, Box const &bbox) {
-	fast_clip(polygon.outer(), bbox);
+	Ring result;
+	fast_clip(polygon.outer(), bbox, result);
 	if (polygon.outer().empty()) {
 		polygon.inners().resize(0);
 		return;
 	}
 	for (auto &inner: polygon.inners()) {
-		fast_clip(inner, bbox);
+		fast_clip(inner, bbox, result);
 	}
 	polygon.inners().erase(std::remove_if(
 		polygon.inners().begin(), polygon.inners().end(), 

src/osm_lua_processing.cpp

@@ -512,21 +512,43 @@ void reverse_project(DegPoint& p) {
     geom::set<1>(p, latp2lat(geom::get<1>(p)));
 }
 
+template <typename DstRing, typename SrcRing>
+void projectRing(DstRing& dst, const SrcRing& src) {
+	dst.resize(src.size());
+	for (std::size_t i = 0; i < src.size(); ++i) {
+		geom::set<0>(dst[i], geom::get<0>(src[i]));
+		geom::set<1>(dst[i], latp2lat(geom::get<1>(src[i])));
+	}
+}
+
+#if BOOST_VERSION >= 106700
+double OsmLuaProcessing::projectedPolygonArea(const Polygon &p) {
+	areaPolygonCache.inners().resize(p.inners().size());
+	projectRing(areaPolygonCache.outer(), p.outer());
+	for (std::size_t i = 0; i < p.inners().size(); ++i) {
+		projectRing(areaPolygonCache.inners()[i], p.inners()[i]);
+	}
+
+	geom::strategy::area::spherical<> sph_strategy(RadiusMeter);
+	return geom::area(areaPolygonCache, sph_strategy);
+}
+#endif
+
 // Returns area
 double OsmLuaProcessing::Area() {
 	if (!IsClosed()) return 0;
 
 #if BOOST_VERSION >= 106700
-	geom::strategy::area::spherical<> sph_strategy(RadiusMeter);
 	if (isRelation) {
 		// Boost won't calculate area of a multipolygon, so we just total up the member polygons
-		return multiPolygonArea(multiPolygonCached());
+		double totalArea = 0;
+		const MultiPolygon &mp = multiPolygonCached();
+		for (MultiPolygon::const_iterator it = mp.begin(); it != mp.end(); ++it) {
+			totalArea += projectedPolygonArea(*it);
+		}
+		return totalArea;
 	} else if (isWay) {
-		// Reproject back into lat/lon and then run Boo
-		geom::model::polygon<DegPoint> p;
-		geom::assign(p,polygonCached());
-		geom::for_each_point(p, reverse_project);
-		return geom::area(p, sph_strategy);
+		return projectedPolygonArea(polygonCached());
 	}
 #else
 	if (isRelation) {
@@ -652,10 +674,12 @@ void OsmLuaProcessing::Layer(const string &layerName, bool area) {
 			}
 			else if (isWay) {
 				//Is there a more efficient way to do this?
-				Linestring ls = linestringCached();
+				const Linestring &ls = linestringCached();
 				Polygon p;
+				p.outer().reserve(ls.size());
 				geom::assign_points(p, ls);
-				mp.push_back(p);
+				mp.reserve(1);
+				mp.push_back(std::move(p));
 
 				auto correctionResult = CorrectGeometry(mp);
 				if(correctionResult == CorrectGeometryResult::Invalid) return;
@@ -872,8 +896,10 @@ Point OsmLuaProcessing::calculateCentroid(CentroidAlgorithm algorithm) {
 				geom::centroid(ls, centroid);
 			}
 		} else {
+			const Linestring &ls = linestringCached();
 			Polygon p;
-			geom::assign_points(p, linestringCached());
+			p.outer().reserve(ls.size());
+			geom::assign_points(p, ls);
 
 			if (algorithm == CentroidAlgorithm::Polylabel) {
 				// CONSIDER: pick precision intelligently