[test] improve point sampler test

test/sampling/test_pointsampling.cc

 #include <stdexcept>
 #include <string>
-#include <array>
+#include <vector>
+#include <cmath>
 
 #include <BRepTools.hxx>
 #include <BRep_Builder.hxx>
@@ -8,51 +9,89 @@
 
 #include <frackit/sampling/makeuniformpointsampler.hh>
 #include <frackit/geometry/point.hh>
+#include <frackit/geometry/box.hh>
+#include <frackit/geometry/boundingbox.hh>
 #include <frackit/geometry/cylinder.hh>
 #include <frackit/geometry/hollowcylinder.hh>
+#include <frackit/geometryutilities/getboundingbox.hh>
 #include <frackit/occ/breputilities.hh>
 
-//! test random sampling of points on geometries
-int main()
+//! sample points on a geometry
+//! if a brep file name is provided, the result is written
+//! into a brep file, such that the result can be checked visually
+template<class G>
+void testPointSampling(const G& geometry,
+                       std::size_t numSamples,
+                       const std::string& brepFileName = "")
 {
-    using ctype = double;
-
     using namespace Frackit;
-    using Cylinder = Cylinder<ctype>;
+    using ctype = typename G::ctype;
+    static constexpr int worldDim = G::worldDimension();
 
-    Cylinder cylinder(0.5, 1.0);
-    auto cylPointSampler = makeUniformPointSampler(cylinder);
+    // sample points
+    auto pointSampler = makeUniformPointSampler(geometry);
+    std::vector<Frackit::Point<ctype, worldDim>> points(numSamples);
+    for (unsigned int i = 0; i < numSamples; ++i) points[i] = pointSampler();
 
-    // sample 500 points
-    std::array<Frackit::Point<ctype, 3>, 500> points;
-    for (unsigned int i = 0; i < 500; ++i)
-        points[i] = cylPointSampler();
+    // check if all of the points are contained in the geometry
+    if ( std::any_of(points.begin(), points.end(),
+                     [&geometry] (const auto& p)
+                     { return !geometry.contains(p); }) )
+        throw std::runtime_error("Point not contained in geometry");
 
-    // create a single compound and write to .brep file
-    // build a single compound shape
-    BRep_Builder b;
-    TopoDS_Compound c;
-    b.MakeCompound(c);
+    // compute bounding box of the point cloud
+    BoundingBox<ctype> bbpc;
     for (const auto& p : points)
-        b.Add(c, OCCUtilities::getShape(p));
-    b.Add(c, OCCUtilities::getShape(cylinder));
+        bbpc += getBoundingBox(p);
 
-    BRepTools::Write(c, "cylinderpoints.brep");
+    // we expect less than 5% deviation to the geometry's bbox
+    const auto bbg = getBoundingBox(geometry);
+    const auto dx = bbg.xMax() - bbg.xMin();
+    const auto dy = bbg.yMax() - bbg.yMin();
+    const auto dz = bbg.zMax() - bbg.zMin();
 
-    // test sampling on a hollow cylinder
-    HollowCylinder<ctype> hollowCylinder(0.4, 0.5, 1.0);
-    cylPointSampler = makeUniformPointSampler(hollowCylinder);
-    for (unsigned int i = 0; i < 500; ++i)
-        points[i] = cylPointSampler();
+    using std::abs;
+    if ( abs(bbpc.xMin() - bbg.xMin()) > 0.05*dx ) throw std::runtime_error("Bounding box deviation > 5%");
+    if ( abs(bbpc.xMax() - bbg.xMax()) > 0.05*dx ) throw std::runtime_error("Bounding box deviation > 5%");
+    if ( abs(bbpc.yMin() - bbg.yMin()) > 0.05*dy ) throw std::runtime_error("Bounding box deviation > 5%");
+    if ( abs(bbpc.yMax() - bbg.yMax()) > 0.05*dy ) throw std::runtime_error("Bounding box deviation > 5%");
+    if ( abs(bbpc.zMin() - bbg.zMin()) > 0.05*dz ) throw std::runtime_error("Bounding box deviation > 5%");
+    if ( abs(bbpc.zMax() - bbg.zMax()) > 0.05*dz ) throw std::runtime_error("Bounding box deviation > 5%");
 
-    // create a single compound and write to .brep file
-    // build a single compound shape
-    TopoDS_Compound c2;
-    b.MakeCompound(c2);
-    for (const auto& p : points)
-        b.Add(c2, OCCUtilities::getShape(p));
-    b.Add(c2, OCCUtilities::getShape(cylinder));
-    BRepTools::Write(c2, "cylinderpoints_hollow.brep");
+    std::cout << "Successfully tested geometry: " << geometry.name() << std::endl;
+    std::cout << "\t- Deviation in x-direction: " << abs(dx - bbpc.xMax() + bbpc.xMin())/dx*100.0 << "%" << std::endl;
+    std::cout << "\t- Deviation in y-direction: " << abs(dy - bbpc.yMax() + bbpc.yMin())/dy*100.0 << "%" << std::endl;
+    std::cout << "\t- Deviation in z-direction: " << abs(dz - bbpc.zMax() + bbpc.zMin())/dz*100.0 << "%" << std::endl;
+
+    if (brepFileName != "")
+    {
+        // create a single compound and write to .brep file
+        BRep_Builder b;
+        TopoDS_Compound c;
+        b.MakeCompound(c);
+        for (const auto& p : points)
+            b.Add(c, OCCUtilities::getShape(p));
+        b.Add(c, OCCUtilities::getShape(geometry));
+
+        BRepTools::Write(c, brepFileName.c_str());
+    }
+}
+
+//! test random sampling of points on geometries
+int main()
+{
+    using ctype = double;
+
+    using namespace Frackit;
+    using Box = Box<ctype>;
+    using BBox = BoundingBox<ctype>;
+    using Cylinder = Cylinder<ctype>;
+    using HollowCyl = HollowCylinder<ctype>;
+
+    testPointSampling(Box(0.,0.,0.,1.,1.,1.), 1000, "pointsamples_box.brep");
+    testPointSampling(BBox(0.,0.,0.,1.,1.,1.), 1000, "pointsamples_bbox.brep");
+    testPointSampling(Cylinder(0.5, 1.0), 1000, "pointsamples_cyl.brep");
+    testPointSampling(HollowCyl(0.4, 0.5, 1.0), 1000, "pointsamples_hollowcyl.brep");
 
     std::cout << "All tests passed" << std::endl;
     return 0;