diff --git a/Octree/include/CGAL/Octree.h b/Octree/include/CGAL/Octree.h index 6bb37663e8c..67deccafd78 100644 --- a/Octree/include/CGAL/Octree.h +++ b/Octree/include/CGAL/Octree.h @@ -155,7 +155,7 @@ namespace CGAL { m_unit_per_depth.push_back(1 << (m_max_depth_reached - i)); } - void refine() { + void _refine() { } @@ -178,6 +178,117 @@ namespace CGAL { } + private: // functions : + + Point compute_barycenter_position(Node &node) const { + + // Determine the side length of this node + FT size = m_side_per_depth[node.depth()]; + + // Determine the location this node should be split + // TODO: I think Point_3 has a [] operator, so using an array here might not be necessary! + + FT bary[3]; + for (int i = 0; i < 3; i++) + bary[i] = node.location()[i] * size + (size / 2.0) + m_bbox_min[i]; + + // Convert that location into a point + return {bary[0], bary[1], bary[2]}; + } + + void refine_recurse(Node &node, size_t dist_to_max_depth, size_t max_pts_num) { + + // Check if the depth limit is reached, or if the node isn't filled + if (dist_to_max_depth == 0 || node.num_points() <= max_pts_num) { + + // If this node is the deepest in the tree, record its depth + if (m_max_depth_reached < node.depth()) m_max_depth_reached = node.depth(); + + // Don't split this node + return; + } + + // Create child nodes + node.split(); + + // Distribute this nodes points among its children + reassign_points(node); + + // Repeat this process for all children (recursive) + for (int child_id = 0; child_id < 8; child_id++) { + refine_recurse(node[child_id], dist_to_max_depth - 1, max_pts_num); + } + } + + void _refine_recurse(Node &node, size_t dist_to_max_depth, size_t max_pts_num) { + + // Check if the depth limit is reached, or if the node isn't filled + if (dist_to_max_depth == 0 || node.num_points() <= max_pts_num) { + + // If this node is the deepest in the tree, record its depth + if (m_max_depth_reached < node.depth()) m_max_depth_reached = node.depth(); + + // Don't split this node + return; + } + + // Create child nodes + node.split(); + + // Distribute this nodes points among its children + _reassign_points(node); + + // Repeat this process for all children (recursive) + for (int child_id = 0; child_id < 8; child_id++) { + _refine_recurse(node[child_id], dist_to_max_depth - 1, max_pts_num); + } + } + + void reassign_points(Node &node) { + + // Find the position of this node's split + Point barycenter = compute_barycenter_position(node); + + // Check each point contained by this node + for (const InputIterator &pwn_it : node.points()) { + const Point &point = get(m_points_map, *pwn_it); + + // Determine which octant a point falls in + // TODO: Could this use std::bitset? + int is_right = barycenter[0] < point[0]; + int is_up = barycenter[1] < point[1]; + int is_front = barycenter[2] < point[2]; + + // Check if a point is very close to the edge + bool equal_right = std::abs(barycenter[0] - point[0]) < 1e-6; + bool equal_up = std::abs(barycenter[1] - point[1]) < 1e-6; + bool equal_front = std::abs(barycenter[2] - point[2]) < 1e-6; + + // Generate a 3-bit code representing a point's octant + int child_id = (is_front << 2) | (is_up << 1) | is_right; + + // Get the child node using that code, and add the point + node[child_id].add_point(pwn_it); + + // Edge cases get special treatment to prevent extremely deep trees + + if (equal_right) { + int sym_child_id = (is_front << 2) | (is_up << 1) | (!is_right); + node[sym_child_id].add_point(pwn_it); + } + + if (equal_up) { + int sym_child_id = (is_front << 2) | (!is_up << 1) | is_right; + node[sym_child_id].add_point(pwn_it); + } + + if (equal_front) { + int sym_child_id = (!is_front << 2) | (is_up << 1) | (!is_right); + node[sym_child_id].add_point(pwn_it); + } + } + } + std::array partition_around_point(Range_iterator begin, Range_iterator end, Point point) { auto partitions = std::array(); @@ -230,104 +341,16 @@ namespace CGAL { return partitions; } - void distribute_points(Node &node, Point center) { + void _reassign_points(Node &node) { - if (node.is_leaf()) return; + Point center = compute_barycenter_position(node); auto partitions = partition_around_point(node._m_points_begin, node._m_points_end, center); + for (int i = 0; i < 8; ++i) { - for (auto it = partitions[i]; it != partitions[i+1]; ++it) - std::cout << *it << ", "; - std::cout << std::endl; - } - } - - - private: // functions : - - Point compute_barycenter_position(Node &node) const { - - // Determine the side length of this node - FT size = m_side_per_depth[node.depth()]; - - // Determine the location this node should be split - // TODO: I think Point_3 has a [] operator, so using an array here might not be necessary! - - FT bary[3]; - for (int i = 0; i < 3; i++) - bary[i] = node.location()[i] * size + (size / 2.0) + m_bbox_min[i]; - - // Convert that location into a point - return {bary[0], bary[1], bary[2]}; - } - - void refine_recurse(Node &node, size_t dist_to_max_depth, size_t max_pts_num) { - - // Check if the depth limit is reached, or if the node isn't filled - if (dist_to_max_depth == 0 || node.num_points() <= max_pts_num) { - - // If this node is the deepest in the tree, record its depth - if (m_max_depth_reached < node.depth()) m_max_depth_reached = node.depth(); - - // Don't split this node - return; - } - - // Create child nodes - node.split(); - - // Distribute this nodes points among its children - reassign_points(node); - - // Repeat this process for all children (recursive) - for (int child_id = 0; child_id < 8; child_id++) { - refine_recurse(node[child_id], dist_to_max_depth - 1, max_pts_num); - } - } - - void reassign_points(Node &node) { - - // Find the position of this node's split - Point barycenter = compute_barycenter_position(node); - - // Check each point contained by this node - for (const InputIterator &pwn_it : node.points()) { - const Point &point = get(m_points_map, *pwn_it); - - // Determine which octant a point falls in - // TODO: Could this use std::bitset? - int is_right = barycenter[0] < point[0]; - int is_up = barycenter[1] < point[1]; - int is_front = barycenter[2] < point[2]; - - // Check if a point is very close to the edge - bool equal_right = std::abs(barycenter[0] - point[0]) < 1e-6; - bool equal_up = std::abs(barycenter[1] - point[1]) < 1e-6; - bool equal_front = std::abs(barycenter[2] - point[2]) < 1e-6; - - // Generate a 3-bit code representing a point's octant - int child_id = (is_front << 2) | (is_up << 1) | is_right; - - // Get the child node using that code, and add the point - node[child_id].add_point(pwn_it); - - // Edge cases get special treatment to prevent extremely deep trees - - if (equal_right) { - int sym_child_id = (is_front << 2) | (is_up << 1) | (!is_right); - node[sym_child_id].add_point(pwn_it); - } - - if (equal_up) { - int sym_child_id = (is_front << 2) | (!is_up << 1) | is_right; - node[sym_child_id].add_point(pwn_it); - } - - if (equal_front) { - int sym_child_id = (!is_front << 2) | (is_up << 1) | (!is_right); - node[sym_child_id].add_point(pwn_it); - } + node[i]._m_points_begin = partitions[i]; + node[i]._m_points_end = partitions[i + 1]; } } diff --git a/Octree/include/CGAL/Octree/Octree_node.h b/Octree/include/CGAL/Octree/Octree_node.h index f3f586d9a18..0c7695a249f 100644 --- a/Octree/include/CGAL/Octree/Octree_node.h +++ b/Octree/include/CGAL/Octree/Octree_node.h @@ -141,6 +141,10 @@ namespace CGAL { size_t num_points() const { return m_points.size(); } + size_t _num_points() const { + return std::distance(_m_points_begin, _m_points_end); + } + bool is_empty() const { return (m_points.size() == 0); } IntPoint &location() { return m_location; }