diff --git a/src/spatialjoin/GeometryCache.cpp b/src/spatialjoin/GeometryCache.cpp index 6b6d3af..6f3cc27 100644 --- a/src/spatialjoin/GeometryCache.cpp +++ b/src/spatialjoin/GeometryCache.cpp @@ -183,6 +183,9 @@ sj::Line sj::GeometryCache::getFromDisk(size_t off, // OBB readPoly(_geomsFReads[tid], ret.obb); + // convex hull + readPoly(_geomsFReads[tid], ret.convexHull); + return ret; } @@ -244,6 +247,9 @@ sj::Area sj::GeometryCache::getFromDisk(size_t off, // OBB readPoly(_geomsFReads[tid], ret.obb); + // convex hull + readPoly(_geomsFReads[tid], ret.convexHull); + // simplified inner readPoly(_geomsFReads[tid], ret.inner); @@ -414,6 +420,9 @@ size_t sj::GeometryCache::writeTo(const sj::Line& val, // OBB ret += writePoly(val.obb, str); + // convex hull + ret += writePoly(val.convexHull, str); + return ret; } @@ -476,6 +485,9 @@ size_t sj::GeometryCache::writeTo(const sj::Area& val, // OBB ret += writePoly(val.obb, str); + // convex hull + ret += writePoly(val.convexHull, str); + // innerGeom ret += writePoly(val.inner, str); diff --git a/src/spatialjoin/GeometryCache.h b/src/spatialjoin/GeometryCache.h index eb2fa17..9bb094f 100644 --- a/src/spatialjoin/GeometryCache.h +++ b/src/spatialjoin/GeometryCache.h @@ -50,6 +50,9 @@ struct Area { // cutouts std::map cutouts; + // Convex hull + util::geo::I32XSortedPolygon convexHull; + // OBB util::geo::I32XSortedPolygon obb; @@ -103,6 +106,9 @@ struct Line { // cutouts std::map cutouts; + // Convex hull + util::geo::I32XSortedPolygon convexHull; + // OBB util::geo::I32XSortedPolygon obb; }; diff --git a/src/spatialjoin/SpatialJoinMain.cpp b/src/spatialjoin/SpatialJoinMain.cpp index 8902007..59e7ac5 100755 --- a/src/spatialjoin/SpatialJoinMain.cpp +++ b/src/spatialjoin/SpatialJoinMain.cpp @@ -74,6 +74,8 @@ void printHelp(int argc, char** argv) { << "disable cutouts\n" << std::setw(41) << " --no-diag-box" << "disable diagonal bounding-box based pre-filter\n" + << std::setw(41) << " --no-convex-hulls" + << "disable convex hulls\n" << std::setw(41) << " --no-fast-sweep-skip" << "disable fast sweep skip using binary search\n" << std::setw(41) << " --use-inner-outer" @@ -118,6 +120,7 @@ int main(int argc, char** argv) { bool useOBB = true; bool useCutouts = true; bool useDiagBox = true; + bool useConvexHulls = true; bool useFastSweepSkip = true; bool useInnerOuter = false; @@ -170,6 +173,8 @@ int main(int argc, char** argv) { useCutouts = false; } else if (cur == "--no-diag-box") { useDiagBox = false; + } else if (cur == "--no-convex-hulls") { + useConvexHulls = false; } else if (cur == "--no-fast-sweep-skip") { useFastSweepSkip = false; } else if (cur == "--use-inner-outer") { @@ -245,7 +250,7 @@ int main(int argc, char** argv) { Sweeper sweeper( {numThreads, numCaches, prefix, intersects, contains, covers, touches, equals, overlaps, crosses, suffix, useBoxIds, useArea, useOBB, - useCutouts, useDiagBox, useFastSweepSkip, useInnerOuter}, + useCutouts, useDiagBox, useConvexHulls, useFastSweepSkip, useInnerOuter}, cache, output); LOGTO(INFO, std::cerr) << "Parsing input geometries..."; diff --git a/src/spatialjoin/Stats.h b/src/spatialjoin/Stats.h index 2c6b656..f0d6673 100644 --- a/src/spatialjoin/Stats.h +++ b/src/spatialjoin/Stats.h @@ -23,6 +23,11 @@ struct Stats { uint64_t timeOBBIsectAreaPoint = 0; uint64_t timeOBBIsectLineLine = 0; + uint64_t timeConvexHullIsectAreaArea = 0; + uint64_t timeConvexHullIsectAreaLine = 0; + uint64_t timeConvexHullIsectAreaPoint = 0; + uint64_t timeConvexHullIsectLineLine = 0; + uint64_t timeFullGeoCheckAreaArea = 0; uint64_t timeFullGeoCheckAreaLine = 0; uint64_t timeFullGeoCheckAreaPoint = 0; @@ -69,7 +74,9 @@ inline std::string Stats::toString() { timeGeoCacheRetrievalPoint + timeWrite + timeBoxIdIsectAreaArea + timeBoxIdIsectAreaLine + timeOBBIsectAreaArea + timeOBBIsectAreaLine + timeOBBIsectAreaPoint + - timeOBBIsectLineLine + timeBoxIdIsectAreaPoint + + timeOBBIsectLineLine + timeConvexHullIsectAreaArea + + timeConvexHullIsectAreaLine + timeConvexHullIsectAreaPoint + + timeConvexHullIsectLineLine + timeBoxIdIsectAreaPoint + timeBoxIdIsectLineLine + timeBoxIdIsectLinePoint + timeFullGeoCheckAreaArea + timeFullGeoCheckAreaLine + timeFullGeoCheckAreaPoint + timeFullGeoCheckLineLine + @@ -129,6 +136,38 @@ inline std::string Stats::toString() { ss << "time for obb intersections LINE/LINE: " << t << " s (" << ((t / sum) * 100.0) << "%)\n"; + t = double(timeConvexHullIsectAreaArea) / 1000000000.0; + ss << "time for convex hull intersections AREA/AREA: " << t << " s (" + << ((t / sum) * 100.0) << "%)\n"; + + t = double(timeConvexHullIsectAreaLine) / 1000000000.0; + ss << "time for convex hull intersections AREA/LINE: " << t << " s (" + << ((t / sum) * 100.0) << "%)\n"; + + t = double(timeConvexHullIsectAreaPoint) / 1000000000.0; + ss << "time for convex hull intersections AREA/POINT: " << t << " s (" + << ((t / sum) * 100.0) << "%)\n"; + + t = double(timeConvexHullIsectLineLine) / 1000000000.0; + ss << "time for convex hull intersections LINE/LINE: " << t << " s (" + << ((t / sum) * 100.0) << "%)\n"; + + t = double(timeConvexHullIsectAreaArea) / 1000000000.0; + ss << "time for convex hull intersections AREA/AREA: " << t << " s (" + << ((t / sum) * 100.0) << "%)\n"; + + t = double(timeConvexHullIsectAreaLine) / 1000000000.0; + ss << "time for convex hull intersections AREA/LINE: " << t << " s (" + << ((t / sum) * 100.0) << "%)\n"; + + t = double(timeConvexHullIsectAreaPoint) / 1000000000.0; + ss << "time for convex hull intersections AREA/POINT: " << t << " s (" + << ((t / sum) * 100.0) << "%)\n"; + + t = double(timeConvexHullIsectLineLine) / 1000000000.0; + ss << "time for convex hull intersections LINE/LINE: " << t << " s (" + << ((t / sum) * 100.0) << "%)\n"; + t = double(timeFullGeoCheckAreaArea) / 1000000000.0; ss << "time for " << fullGeoChecksAreaArea << " full geom checks AREA/AREA: " << t << " s (" << ((t / sum) * 100.0) @@ -280,6 +319,10 @@ inline Stats operator+(const Stats& a, const Stats& b) { a.timeOBBIsectAreaLine + b.timeOBBIsectAreaLine, a.timeOBBIsectAreaPoint + b.timeOBBIsectAreaPoint, a.timeOBBIsectLineLine + b.timeOBBIsectLineLine, + a.timeConvexHullIsectAreaArea + b.timeConvexHullIsectAreaArea, + a.timeConvexHullIsectAreaLine + b.timeConvexHullIsectAreaLine, + a.timeConvexHullIsectAreaPoint + b.timeConvexHullIsectAreaPoint, + a.timeConvexHullIsectLineLine + b.timeConvexHullIsectLineLine, a.timeFullGeoCheckAreaArea + b.timeFullGeoCheckAreaArea, a.timeFullGeoCheckAreaLine + b.timeFullGeoCheckAreaLine, a.timeFullGeoCheckAreaPoint + b.timeFullGeoCheckAreaPoint, diff --git a/src/spatialjoin/Sweeper.cpp b/src/spatialjoin/Sweeper.cpp index 836cf12..7c39941 100644 --- a/src/spatialjoin/Sweeper.cpp +++ b/src/spatialjoin/Sweeper.cpp @@ -46,6 +46,7 @@ using util::geo::webMercToLatLng; const static size_t CUTOUTS_MIN_SIZE = 100; const static size_t OBB_MIN_SIZE = 100; +const static size_t CONVEXHULL_MIN_SIZE = 100; const static double sin45 = 1.0 / sqrt(2); const static double cos45 = 1.0 / sqrt(2); @@ -136,6 +137,7 @@ void Sweeper::add(const I32Polygon& poly, const std::string& gid, size_t subid, WriteBatch& batch) { WriteCand cur; const auto& box = getBoundingBox(poly); + const auto& hull = util::geo::convexHull(poly); I32XSortedPolygon spoly(poly); double areaSize = area(poly); @@ -214,7 +216,6 @@ void Sweeper::add(const I32Polygon& poly, const std::string& gid, size_t subid, } util::geo::I32Polygon obb; - if (_cfg.useOBB && poly.getOuter().size() >= OBB_MIN_SIZE) { obb = util::geo::convexHull( util::geo::pad(util::geo::getOrientedEnvelope(poly), 10)); @@ -223,11 +224,29 @@ void Sweeper::add(const I32Polygon& poly, const std::string& gid, size_t subid, if (obb.getOuter().size() >= poly.getOuter().size()) obb = {}; } + util::geo::I32Polygon convexHull; + if (_cfg.useConvexHull && poly.getOuter().size() >= CONVEXHULL_MIN_SIZE) { + convexHull = util::geo::convexHull(poly); + + // drop redundant convex hull + if (convexHull.getOuter().size() >= poly.getOuter().size()) + convexHull = {}; + if (!obb.getOuter().empty() && convexHull.getOuter().size() == obb.getOuter().size() && outerArea(convexHull) == outerArea(obb)) + convexHull = {}; + if (!outer.empty() && convexHull.getOuter().size() == outer.getOuter().rawRing().size()) + convexHull = {}; + } + + // if (subid > 0) { + // std::unique_lock lock(_multiAddMtx); + // multiAdd(gid, box.getLowerLeft().getX(), box.getUpperRight().getX()); + // } + std::stringstream str; GeometryCache::writeTo( {spoly, box, gid, subid, areaSize, _cfg.useArea ? outerAreaSize : 0, - boxIds, cutouts, obb, inner, innerBox, innerOuterAreaSize, outer, - outerBox, outerOuterAreaSize}, + boxIds, cutouts, obb, convexHull, inner, innerBox, innerOuterAreaSize, + outer, outerBox, outerOuterAreaSize}, str); ; @@ -327,13 +346,23 @@ void Sweeper::add(const I32Line& line, const std::string& gid, size_t subid, if (obb.getOuter().size() >= line.size()) obb = {}; } + util::geo::I32Polygon convexHull; + if (_cfg.useConvexHull && line.size() >= CONVEXHULL_MIN_SIZE) { + convexHull = util::geo::convexHull(line); + + // drop redundant hull + if (convexHull.getOuter().size() >= line.size()) convexHull = {}; + if (!obb.getOuter().empty() && convexHull.getOuter().size() == obb.getOuter().size() && outerArea(convexHull) == outerArea(obb)) + convexHull = {}; + } + if (!_cfg.useFastSweepSkip) { sline.setMaxSegLen(std::numeric_limits::max()); } std::stringstream str; GeometryCache::writeTo( - {sline, box, gid, subid, len, boxIds, cutouts, obb}, str); + {sline, box, gid, subid, len, boxIds, cutouts, obb, convexHull}, str); cur.raw = str.str(); cur.boxvalIn = {0, // placeholder, will be overwritten later on @@ -1003,6 +1032,7 @@ sj::Area Sweeper::areaFromSimpleArea(const SimpleArea* sa) const { {}, {}, {}, + {}, 0, {}, {}, @@ -1054,13 +1084,19 @@ std::tuple Sweeper::check(const Area* a, } } - if (_cfg.useOBB) { - if (!a->obb.empty() && !b->obb.empty()) { - auto ts = TIME(); - auto r = util::geo::intersectsContainsCovers(a->obb, b->obb); - _stats[t].timeOBBIsectAreaArea += TOOK(ts); - if (!std::get<0>(r)) return {0, 0, 0, 0, 0}; - } + if (_cfg.useOBB && !a->obb.empty() && !b->obb.empty()) { + auto ts = TIME(); + auto r = util::geo::intersectsContainsCovers(a->obb, b->obb); + _stats[t].timeOBBIsectAreaArea += TOOK(ts); + if (!std::get<0>(r)) return {0, 0, 0, 0, 0}; + } + + if (_cfg.useConvexHull && !a->convexHull.empty() && !b->convexHull.empty()) { + auto ts = TIME(); + auto r = util::geo::intersectsContainsCovers(a->convexHull, a->box, 0, + b->convexHull, b->box, 0); + _stats[t].timeConvexHullIsectAreaArea += TOOK(ts); + if (!std::get<0>(r)) return {0, 0, 0, 0, 0}; } if (_cfg.useInnerOuter && !a->outer.empty() && !b->outer.empty()) { @@ -1160,6 +1196,14 @@ std::tuple Sweeper::check(const Line* a, } } + if (_cfg.useConvexHull && !a->convexHull.empty() && !b->convexHull.empty()) { + auto ts = TIME(); + auto r = intersectsContainsCovers(a->convexHull, a->box, 0, b->convexHull, + b->box, 0); + _stats[t].timeConvexHullIsectAreaLine += TOOK(ts); + if (!std::get<0>(r)) return {0, 0, 0, 0, 0}; + } + if (_cfg.useInnerOuter && !b->outer.empty()) { auto ts = TIME(); auto r = util::geo::intersectsContainsCovers(a->geom, a->box, b->outer, @@ -1229,13 +1273,19 @@ std::tuple Sweeper::check(const Line* a, } } - if (_cfg.useOBB) { - if (!a->obb.empty() && !b->obb.empty()) { - auto ts = TIME(); - auto r = util::geo::intersectsContainsCovers(a->obb, b->obb); - _stats[t].timeOBBIsectLineLine += TOOK(ts); - if (!std::get<0>(r)) return {0, 0, 0, 0, 0}; - } + if (_cfg.useOBB && !a->obb.empty() && !b->obb.empty()) { + auto ts = TIME(); + auto r = util::geo::intersectsContainsCovers(a->obb, b->obb); + _stats[t].timeOBBIsectLineLine += TOOK(ts); + if (!std::get<0>(r)) return {0, 0, 0, 0, 0}; + } + + if (_cfg.useConvexHull && !a->convexHull.empty() && !b->convexHull.empty()) { + auto ts = TIME(); + auto r = intersectsContainsCovers(a->convexHull, a->box, 0, b->convexHull, + b->box, 0); + _stats[t].timeConvexHullIsectLineLine += TOOK(ts); + if (!std::get<0>(r)) return {0, 0, 0, 0, 0}; } auto ts = TIME(); @@ -1288,6 +1338,14 @@ std::tuple Sweeper::check(const SimpleLine* a, if (!std::get<0>(r)) return {0, 0, 0, 0, 0}; } + if (_cfg.useConvexHull && !b->convexHull.empty()) { + auto ts = TIME(); + auto r = intersectsContainsCovers( + I32XSortedLine(LineSegment(a->a, a->b)), b->convexHull); + _stats[t].timeConvexHullIsectLineLine += TOOK(ts); + if (!std::get<0>(r)) return {0, 0, 0, 0, 0}; + } + if (_cfg.useInnerOuter && !b->outer.empty()) { auto ts = TIME(); auto r = util::geo::intersectsContainsCovers( @@ -1476,6 +1534,13 @@ std::pair Sweeper::check(const I32Point& a, const Area* b, if (!std::get<1>(r)) return {0, 0}; } + if (_cfg.useConvexHull && !b->convexHull.empty()) { + auto ts = TIME(); + auto r = containsCovers(a, b->convexHull); + _stats[t].timeConvexHullIsectAreaPoint += TOOK(ts); + if (!std::get<1>(r)) return {0, 0}; + } + if (_cfg.useInnerOuter && !b->outer.empty()) { auto ts = TIME(); auto r = containsCovers(a, b->outer); diff --git a/src/spatialjoin/Sweeper.h b/src/spatialjoin/Sweeper.h index 065ff8a..8ff5788 100644 --- a/src/spatialjoin/Sweeper.h +++ b/src/spatialjoin/Sweeper.h @@ -109,6 +109,7 @@ struct SweeperCfg { bool useDiagBox; bool useFastSweepSkip; bool useInnerOuter; + bool useConvexHull; }; // buffer size _must_ be multiples of sizeof(BoxVal) diff --git a/src/spatialjoin/tests/TestMain.cpp b/src/spatialjoin/tests/TestMain.cpp index 95b13ff..5769c83 100644 --- a/src/spatialjoin/tests/TestMain.cpp +++ b/src/spatialjoin/tests/TestMain.cpp @@ -67,59 +67,66 @@ int main(int, char**) { 1, 1, "$", " intersects ", " contains ", " covers ", " touches ", " equals ", " overlaps ", " crosses ", "$\n", false, false, false, false, false, - false, false}; + false, false, false}; sj::SweeperCfg all{1, 1, "$", " intersects ", " contains ", " covers ", " touches ", " equals ", " overlaps ", " crosses ", "$\n", true, true, true, true, true, true, - true}; + true, true}; sj::SweeperCfg noSurfaceArea{ 1, 1, "$", " intersects ", " contains ", " covers ", " touches ", " equals ", " overlaps ", " crosses ", "$\n", true, false, true, true, true, - true, true}; + true, true, true}; sj::SweeperCfg noBoxIds{ 1, 1, "$", " intersects ", " contains ", " covers ", " touches ", " equals ", " overlaps ", " crosses ", "$\n", false, true, true, true, true, - true, true}; + true, true, true}; sj::SweeperCfg noObb{1, 1, "$", " intersects ", " contains ", " covers ", " touches ", " equals ", " overlaps ", " crosses ", "$\n", true, true, false, true, true, true, - true}; + true, true}; sj::SweeperCfg noCutouts{ 1, 1, "$", " intersects ", " contains ", " covers ", " touches ", " equals ", " overlaps ", " crosses ", "$\n", true, true, true, false, true, - true, true}; + true, true, true}; sj::SweeperCfg noDiagBox{ 1, 1, "$", " intersects ", " contains ", " covers ", " touches ", " equals ", " overlaps ", " crosses ", "$\n", true, true, true, true, false, - true, true}; + true, true, true}; sj::SweeperCfg noFastSweep{ 1, 1, "$", " intersects ", " contains ", " covers ", " touches ", " equals ", " overlaps ", " crosses ", "$\n", true, true, true, true, true, - false, true}; + false, true, true}; sj::SweeperCfg noInnerOuter{ 1, 1, "$", " intersects ", " contains ", " covers ", " touches ", " equals ", " overlaps ", " crosses ", "$\n", true, true, true, true, true, - true, false}; + true, false, true}; + + sj::SweeperCfg noConvexHulls{ + 1, 1, "$", " intersects ", " contains ", " covers ", + " touches ", " equals ", " overlaps ", " crosses ", "$\n", + true, true, true, true, true, + true, true, false}; std::vector cfgs{baseline, all, noSurfaceArea, noBoxIds, noObb, noCutouts, - noDiagBox, noFastSweep, noInnerOuter}; + noDiagBox, noFastSweep, noInnerOuter, + noConvexHulls}; for (auto cfg : cfgs) { {