Convex hull

src/spatialjoin/GeometryCache.cpp

@@ -144,6 +144,9 @@ sj::Line sj::GeometryCache<sj::Line>::getFromDisk(size_t off,
     ret.cutouts[boxid] = cutout;
   }
 
+  // convex hull
+  readPoly(_geomsFReads[tid], ret.convexHull);
+
   // OBB
   readPoly(_geomsFReads[tid], ret.obb);
 
@@ -206,6 +209,9 @@ sj::Area sj::GeometryCache<sj::Area>::getFromDisk(size_t off,
     ret.cutouts[boxid] = cutout;
   }
 
+  // convex hull
+  readPoly(_geomsFReads[tid], ret.convexHull);
+
   // OBB
   readPoly(_geomsFReads[tid], ret.obb);
 
@@ -313,6 +319,9 @@ size_t sj::GeometryCache<sj::Line>::add(const sj::Line& val) {
     _geomsOffset += sizeof(int32_t) + sizeof(int32_t);
   }
 
+  // convex hull
+  writePoly(val.convexHull);
+
   // OBB
   writePoly(val.obb);
 
@@ -375,6 +384,9 @@ size_t sj::GeometryCache<sj::Area>::add(const sj::Area& val) {
     _geomsOffset += sizeof(int32_t) + sizeof(int32_t);
   }
 
+  // convex hull
+  writePoly(val.convexHull);
+
   // OBB
   writePoly(val.obb);
 

src/spatialjoin/GeometryCache.h

@@ -49,6 +49,9 @@ struct Area {
   // cutouts
   std::map<int32_t, size_t> cutouts;
 
+  // Convex hull
+  util::geo::I32XSortedPolygon convexHull;
+
   // OBB
   util::geo::I32XSortedPolygon obb;
 
@@ -90,6 +93,9 @@ struct Line {
   // cutouts
   std::map<int32_t, size_t> cutouts;
 
+  // Convex hull
+  util::geo::I32XSortedPolygon convexHull;
+
   // OBB
   util::geo::I32XSortedPolygon obb;
 };

src/spatialjoin/SpatialJoinMain.cpp

@@ -72,6 +72,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"
@@ -107,6 +109,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;
 
@@ -153,6 +156,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") {
@@ -219,7 +224,7 @@ int main(int argc, char** argv) {
   Sweeper sweeper(
       {NUM_THREADS, prefix, intersects, contains, covers, touches, equals,
        overlaps, crosses, suffix, useBoxIds, useArea, useOBB, useCutouts,
-       useDiagBox, useFastSweepSkip, useInnerOuter},
+       useDiagBox, useConvexHulls, useFastSweepSkip, useInnerOuter},
       useCache, cache, output);
 
   if (!useCache) {

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;
@@ -64,7 +69,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 +
@@ -124,6 +131,22 @@ 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(timeFullGeoCheckAreaArea) / 1000000000.0;
   ss << "time for " << fullGeoChecksAreaArea
      << " full geom checks AREA/AREA: " << t << " s (" << ((t / sum) * 100.0)

src/spatialjoin/Sweeper.cpp

@@ -133,8 +133,8 @@ void Sweeper::add(const I32Polygon& poly, const std::string& gid) {
 void Sweeper::add(const I32Polygon& poly, const std::string& gid,
                   size_t subid) {
   const auto& box = getBoundingBox(poly);
+  const auto& hull = util::geo::convexHull(poly);
   I32XSortedPolygon spoly(poly);
-
   double areaSize = area(poly);
   double outerAreaSize = outerArea(poly);
   BoxIdList boxIds;
@@ -163,6 +163,11 @@ void Sweeper::add(const I32Polygon& poly, const std::string& gid,
     outer = sj::innerouter::simplifiedPoly<Mode::OUTER>(poly, 1 / (3.14 * 20));
   }
 
+  util::geo::I32Polygon convexHull;
+  if (_cfg.useConvexHull) {
+    convexHull = util::geo::convexHull(poly);
+  }
+
   util::geo::I32Polygon obb;
   obb = util::geo::convexHull(util::geo::getOrientedEnvelope(poly));
 
@@ -178,7 +183,7 @@ void Sweeper::add(const I32Polygon& poly, const std::string& gid,
 
   size_t id = _areaCache.add({spoly, box, gid, subid, areaSize,
                               _cfg.useArea ? outerAreaSize : 0, boxIds, cutouts,
-                              obb, inner, outer});
+                              obb, inner, outer, convexHull});
 
   diskAdd({id, box.getLowerLeft().getY(), box.getUpperRight().getY(),
            box.getLowerLeft().getX(), false, POLYGON, areaSize, box45});
@@ -210,6 +215,11 @@ void Sweeper::add(const I32Line& line, const std::string& gid, size_t subid) {
 
   double len = util::geo::len(line);
 
+  util::geo::I32Polygon convexHull;
+  if (_cfg.useConvexHull) {
+    convexHull = util::geo::convexHull(line);
+  }
+
   util::geo::I32Polygon obb;
   obb = util::geo::convexHull(util::geo::getOrientedEnvelope(line));
 
@@ -241,7 +251,7 @@ void Sweeper::add(const I32Line& line, const std::string& gid, size_t subid) {
     }
 
     size_t id =
-        _lineCache.add({sline, box, gid, subid, len, boxIds, cutouts, obb});
+        _lineCache.add({sline, box, gid, subid, len, boxIds, cutouts, obb, convexHull});
 
     diskAdd({id, box.getLowerLeft().getY(), box.getUpperRight().getY(),
              box.getLowerLeft().getX(), false, LINE, len, box45});
@@ -665,6 +675,15 @@ std::tuple<bool, bool, bool, bool, bool> Sweeper::check(const Area* a,
     }
   }
 
+  if (_cfg.useConvexHull) {
+    auto ts = TIME();
+    auto r = util::geo::intersectsContainsCovers(
+        a->convexHull, a->box, a->outerArea,
+        b->convexHull, b->box, b->outerArea);
+    _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()) {
     auto ts = TIME();
     auto r = util::geo::intersectsContainsCovers(
@@ -759,6 +778,14 @@ std::tuple<bool, bool, bool, bool, bool> Sweeper::check(const Line* a,
     }
   }
 
+  if (_cfg.useConvexHull) {
+    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 =
@@ -839,6 +866,14 @@ std::tuple<bool, bool, bool, bool, bool> Sweeper::check(const Line* a,
     }
   }
 
+  if (_cfg.useConvexHull) {
+    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();
   auto res = intersectsCovers(a->geom, b->geom, a->box, b->box);
   _stats[t].timeFullGeoCheckLineLine += TOOK(ts);
@@ -1081,6 +1116,13 @@ std::pair<bool, bool> Sweeper::check(const I32Point& a, const Area* b,
     if (!std::get<1>(r)) return {0, 0};
   }
 
+  if (_cfg.useConvexHull) {
+    auto ts = TIME();
+    auto r = containsCovers(a, b->convexHull);
+    _stats[t].timeConvexHullIsectAreaPoint += TOOK(ts);
+    if (!std::get<0>(r)) return {0, 0};
+  }
+
   if (_cfg.useInnerOuter && !b->outer.empty()) {
     auto ts = TIME();
     auto r = containsCovers(a, b->outer);

src/spatialjoin/Sweeper.h

@@ -86,6 +86,7 @@ struct SweeperCfg {
   bool useDiagBox;
   bool useFastSweepSkip;
   bool useInnerOuter;
+  bool useConvexHull;
 };
 
 // buffer size _must_ be multiples of sizeof(BoxVal)

src/spatialjoin/tests/TestMain.cpp

@@ -17,7 +17,7 @@ using sj::Sweeper;
 std::string fullRun(const std::string& file) {
   Sweeper sweeper({1, "$", " intersects ", " contains ", " covers ",
                    " touches ", " equals ", " overlaps ", " crosses ", "$\n",
-                   true, true, true, true, true, true, true},
+                   true, true, true, true, true, true, true, true},
                   false, ".", ".resTmp");
 
   size_t gid = 0;