Commit f7b25118 authored by Deepak Rangarajan's avatar Deepak Rangarajan Committed by Jordan Musser
Browse files

try to make hull generation deterministic

parent 78b2a9ea
Loading
Loading
Loading
Loading
Loading
+22 −9
Original line number Diff line number Diff line
@@ -10,7 +10,6 @@
#include <CGAL/boost/graph/graph_traits_Polyhedron_3.h>
#include <CGAL/convex_hull_3.h>
#include <CGAL/point_generators_3.h>
#include <CGAL/Random.h>

namespace {

@@ -19,6 +18,7 @@ typedef typename HalfedgeDS::Vertex Vertex;
typedef typename Vertex::Point Point;
typedef CGAL::Side_of_triangle_mesh<cgal_helper::Polyhedron, cgal_helper::CK>
    Point_inside;
typedef CGAL::Creator_uniform_3<double, cgal_helper::CK::Point_3> Creator;

// A modifier creating a polyhedron with the incremental builder.
template <class HDS> class PolyhedronBuilder : public CGAL::Modifier_base<HDS> {
@@ -60,6 +60,23 @@ public:
  }
};

void points_on_a_sphere(std::list<cgal_helper::CK::Point_3> &points, int N,
                        double radius) {
  points.clear();

  double phi = CGAL_PI * (3.0 - std::sqrt(5.0));
  for (int i = 0; i < N; ++i) {
    double y = 1.0 - (i / double(N - 1)) * 2;
    double r = std::sqrt(1.0 - y * y);
    double theta = phi * i;

    double x = std::cos(theta) * r;
    double z = std::sin(theta) * r;
    cgal_helper::CK::Point_3 pt(radius * x, radius * y, radius * z);
    points.push_back(pt);
  }
}

} // namespace

namespace cgal_helper {
@@ -107,20 +124,16 @@ create_polyhedron_from_hull(const std::tuple<double, double, double> &cube_size,
  // This is currently chosen somewhat arbitrarily
  // based on manual testing and CGAL examples
  // TODO: Is this sufficient?
  const int NUM_SAMPLING_PTS = 250;
  const int NUM_SAMPLING_PTS = 343;

  auto p = std::make_shared<Polyhedron>();

  // Use a fixed seed random function to ensure reproducibility
  auto rand = CGAL::Random(0);

  std::list<cgal_helper::CK::Point_3> points;
  CGAL::Random_points_on_sphere_3<cgal_helper::CK::Point_3> s(sphere_radius, rand);
  std::copy_n(s, NUM_SAMPLING_PTS, std::back_inserter(points));
  points_on_a_sphere(points, NUM_SAMPLING_PTS, sphere_radius);

  CGAL::Random_points_in_cube_3<cgal_helper::CK::Point_3> c(1.0, rand);
  std::list<cgal_helper::CK::Point_3> points_c;
  std::copy_n(c, NUM_SAMPLING_PTS, std::back_inserter(points_c));
  CGAL::points_on_cube_grid_3(1.0, (std::size_t)(NUM_SAMPLING_PTS),
                              std::back_inserter(points_c), Creator());

  auto [sa, sb, sc] = cube_size;

+48 −19
Original line number Diff line number Diff line
@@ -8,8 +8,10 @@ namespace {
TEST_CASE("hull of a sphere and translated cube", "[Levelset Hull]") {
  csg::Hull my_hull;

  double XX = 1.0, YY = 2.0, ZZ = 3.0, Cx = 5.0, Cy = 5.0, Cz = 6.0, R = 1.0;
  double XX = 1.0, YY = 2.0, ZZ = 3.0, Cx = 5.0, Cy = 5.0, Cz = 6.0;

  SECTION("Unit radius"){
     double R = 1.0;
     my_hull.cube_size = {XX, YY, ZZ};
     my_hull.cube_center = {Cx, Cy, Cz};
     my_hull.sphere_radius = R;
@@ -23,11 +25,38 @@ TEST_CASE("hull of a sphere and translated cube", "[Levelset Hull]") {
       CHECK_FALSE(0 < my_levelset(0, 0, Cz));
       CHECK_FALSE(0 < my_levelset(1.3 * R, 0, 0));
       CHECK_FALSE(0 < my_levelset(Cx + 1.1 * XX / 2, Cy, Cz));
       CHECK_FALSE(0 < my_levelset(1.01*(Cx + 0.5*XX), Cy, Cz));
     }
     SECTION("Inside") {
       CHECK(0 < my_levelset(0, 0, 0));
       CHECK(0 < my_levelset(Cx, Cy, Cz));
       CHECK(0 < my_levelset(1.1 * R, 0, 0));
       CHECK(0 < my_levelset(0.99*(Cx + 0.5*XX), Cy, Cz));
     }
  }
  SECTION("Radius 0.1"){
     double R = 0.1;
     my_hull.cube_size = {XX, YY, ZZ};
     my_hull.cube_center = {Cx, Cy, Cz};
     my_hull.sphere_radius = R;
     my_hull.generate_polyhedron();

     auto my_tree = std::make_shared<csg::Tree>();
     my_tree->top.objs.push_back(my_hull);
     csg::CsgIF my_levelset(my_tree);

     SECTION("Outside") {
       CHECK_FALSE(0 < my_levelset(0, 0, Cz));
       CHECK_FALSE(0 < my_levelset(1.4 * R, 0, 0));
       CHECK_FALSE(0 < my_levelset(Cx + 1.1 * XX / 2, Cy, Cz));
       CHECK_FALSE(0 < my_levelset(1.01*(Cx + 0.5*XX), Cy, Cz));
     }
     SECTION("Inside") {
       CHECK(0 < my_levelset(0, 0, 0));
       CHECK(0 < my_levelset(Cx, Cy, Cz));
       CHECK(0 < my_levelset(1.1 * R, 0, 0));
       CHECK(0 < my_levelset(0.99*(Cx + 0.5*XX), Cy, Cz));
     }
  }
}