Geogram  Version 1.9.1
A programming library of geometric algorithms
mesh_AABB.h
Go to the documentation of this file.
1 /*
2  * Copyright (c) 2000-2022 Inria
3  * All rights reserved.
4  *
5  * Redistribution and use in source and binary forms, with or without
6  * modification, are permitted provided that the following conditions are met:
7  *
8  * * Redistributions of source code must retain the above copyright notice,
9  * this list of conditions and the following disclaimer.
10  * * Redistributions in binary form must reproduce the above copyright notice,
11  * this list of conditions and the following disclaimer in the documentation
12  * and/or other materials provided with the distribution.
13  * * Neither the name of the ALICE Project-Team nor the names of its
14  * contributors may be used to endorse or promote products derived from this
15  * software without specific prior written permission.
16  *
17  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
18  * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
19  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
20  * ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
21  * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
22  * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
23  * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
24  * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
25  * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
26  * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
27  * POSSIBILITY OF SUCH DAMAGE.
28  *
29  * Contact: Bruno Levy
30  *
31  * https://www.inria.fr/fr/bruno-levy
32  *
33  * Inria,
34  * Domaine de Voluceau,
35  * 78150 Le Chesnay - Rocquencourt
36  * FRANCE
37  *
38  */
39 
40 #ifndef GEOGRAM_MESH_MESH_AABB
41 #define GEOGRAM_MESH_MESH_AABB
42 
49 #include <geogram/basic/common.h>
50 #include <geogram/mesh/mesh.h>
51 #include <geogram/basic/geometry.h>
52 
53 namespace GEO {
54 
55 
60  template <class BOX> class AABB {
61 
62  protected:
63 
70  void initialize(
71  index_t nb, std::function<void(BOX&, index_t)> get_bbox
72  ) {
73  nb_ = nb;
74  bboxes_.resize(max_node_index(1, 0, nb) + 1);
75  // +1 because size == max_index + 1 !!!
76  init_bboxes_recursive(1, 0, nb_, get_bbox);
77  }
78 
79 
99  std::function<void(index_t)> action, const BOX& box,
100  index_t node, index_t b, index_t e
101  ) const {
102  geo_debug_assert(e != b);
103 
104  // Prune sub-tree that does not have intersection
105  if(!bboxes_overlap(box, bboxes_[node])) {
106  return;
107  }
108 
109  // Leaf case
110  if(e == b+1) {
111  action(b);
112  return;
113  }
114 
115  // Recursion
116  index_t m = b + (e - b) / 2;
117  index_t node_l = 2 * node;
118  index_t node_r = 2 * node + 1;
119 
120  bbox_intersect_recursive(action, box, node_l, b, m);
121  bbox_intersect_recursive(action, box, node_r, m, e);
122  }
123 
146  std::function<void(index_t,index_t)> action,
147  index_t node1, index_t b1, index_t e1,
148  index_t node2, index_t b2, index_t e2
149  ) const {
150  geo_debug_assert(e1 != b1);
151  geo_debug_assert(e2 != b2);
152 
153  // Since we are intersecting the AABBTree with *itself*,
154  // we can prune half of the cases by skipping the test
155  // whenever node2's facet index interval is greated than
156  // node1's facet index interval.
157  if(e2 <= b1) {
158  return;
159  }
160 
161  // The acceleration is here:
162  if(!bboxes_overlap(bboxes_[node1], bboxes_[node2])) {
163  return;
164  }
165 
166  // Simple case: leaf - leaf intersection.
167  if(b1 + 1 == e1 && b2 + 1 == e2) {
168  action(b1, b2);
169  return;
170  }
171 
172  // If node2 has more elements than node1, then
173  // intersect node2's two children with node1
174  // else
175  // intersect node1's two children with node2
176  if(e2 - b2 > e1 - b1) {
177  index_t m2 = b2 + (e2 - b2) / 2;
178  index_t node2_l = 2 * node2;
179  index_t node2_r = 2 * node2 + 1;
180  self_intersect_recursive(action, node1, b1, e1, node2_l, b2, m2);
181  self_intersect_recursive(action, node1, b1, e1, node2_r, m2, e2);
182  } else {
183  index_t m1 = b1 + (e1 - b1) / 2;
184  index_t node1_l = 2 * node1;
185  index_t node1_r = 2 * node1 + 1;
186  self_intersect_recursive(action, node1_l, b1, m1, node2, b2, e2);
187  self_intersect_recursive(action, node1_r, m1, e1, node2, b2, e2);
188  }
189  }
190 
191 
215  std::function<void(index_t,index_t)> action,
216  index_t node1, index_t b1, index_t e1,
217  const AABB<BOX>* other,
218  index_t node2, index_t b2, index_t e2
219  ) const {
220  geo_debug_assert(e1 != b1);
221  geo_debug_assert(e2 != b2);
222 
223  // The acceleration is here:
224  if(!bboxes_overlap(bboxes_[node1], other->bboxes_[node2])) {
225  return;
226  }
227 
228  // Simple case: leaf - leaf intersection.
229  if(b1 + 1 == e1 && b2 + 1 == e2) {
230  action(b1, b2);
231  return;
232  }
233 
234  // If node2 has more elements than node1, then
235  // intersect node2's two children with node1
236  // else
237  // intersect node1's two children with node2
238  if(e2 - b2 > e1 - b1) {
239  index_t m2 = b2 + (e2 - b2) / 2;
240  index_t node2_l = 2 * node2;
241  index_t node2_r = 2 * node2 + 1;
243  action, node1, b1, e1, other, node2_l, b2, m2
244  );
246  action, node1, b1, e1, other, node2_r, m2, e2
247  );
248  } else {
249  index_t m1 = b1 + (e1 - b1) / 2;
250  index_t node1_l = 2 * node1;
251  index_t node1_r = 2 * node1 + 1;
253  action, node1_l, b1, m1, other, node2, b2, e2
254  );
256  action, node1_r, m1, e1, other, node2, b2, e2
257  );
258  }
259  }
260 
261 
269  static index_t max_node_index(index_t node_index, index_t b, index_t e) {
270  geo_debug_assert(e > b);
271  if(b + 1 == e) {
272  return node_index;
273  }
274  index_t m = b + (e - b) / 2;
275  index_t childl = 2 * node_index;
276  index_t childr = 2 * node_index + 1;
277  return std::max(
278  max_node_index(childl, b, m),
279  max_node_index(childr, m, e)
280  );
281  }
282 
294  index_t node_index,
295  index_t b, index_t e,
296  std::function<void(BOX&, index_t)> get_bbox
297  ) {
298  geo_debug_assert(node_index < bboxes_.size());
299  geo_debug_assert(b != e);
300  if(b + 1 == e) {
301  get_bbox(bboxes_[node_index], b);
302  return;
303  }
304  index_t m = b + (e - b) / 2;
305  index_t childl = 2 * node_index;
306  index_t childr = 2 * node_index + 1;
307  geo_debug_assert(childl < bboxes_.size());
308  geo_debug_assert(childr < bboxes_.size());
309  init_bboxes_recursive(childl, b, m, get_bbox);
310  init_bboxes_recursive(childr, m, e, get_bbox);
311  geo_debug_assert(childl < bboxes_.size());
312  geo_debug_assert(childr < bboxes_.size());
313  bbox_union(bboxes_[node_index], bboxes_[childl], bboxes_[childr]);
314  }
315 
316  protected:
317  index_t nb_;
318  vector<BOX> bboxes_;
319  };
320 
321  typedef AABB<Box2d> AABB2d;
322  typedef AABB<Box3d> AABB3d;
323 
324  /**************************************************************/
325 
330  class GEOGRAM_API MeshAABB2d : public AABB2d {
331  public:
335  MeshAABB2d() : mesh_(nullptr) {
336  }
337 
342  const Mesh* mesh() const {
343  return mesh_;
344  }
345 
346  protected:
347  Mesh* mesh_;
348  };
349 
350  /**************************************************************/
351 
356  class GEOGRAM_API MeshAABB3d : public AABB3d {
357  public:
361  MeshAABB3d() : mesh_(nullptr) {
362  }
363 
368  const Mesh* mesh() const {
369  return mesh_;
370  }
371 
372  protected:
373  Mesh* mesh_;
374  };
375 
376  /**************************************************************/
377 
383  class GEOGRAM_API MeshFacetsAABB : public MeshAABB3d {
384  public:
385 
390  struct Intersection {
391  Intersection() :
393  f(index_t(-1)),
394  i(index_t(-1)), j(index_t(-1)), k(index_t(-1))
395  {
396  }
397  vec3 p;
398  double t;
400  vec3 N;
401  index_t i,j,k;
402  double u,v;
403  };
404 
405 
411 
422  void initialize(Mesh& M, bool reorder = true);
423 
424 
435  MeshFacetsAABB(Mesh& M, bool reorder = true);
436 
445  std::function<void(index_t, index_t)> action
446  ) const {
447  self_intersect_recursive(
448  action,
449  1, 0, mesh_->facets.nb(),
450  1, 0, mesh_->facets.nb()
451  );
452  }
453 
462  const Box& box_in, std::function<void(index_t)> action
463  ) const {
464  bbox_intersect_recursive(
465  action, box_in, 1, 0, mesh_->facets.nb()
466  );
467  }
468 
477  const vec3& p, vec3& nearest_point, double& sq_dist
478  ) const {
479  index_t nearest_facet;
480  get_nearest_facet_hint(p, nearest_facet, nearest_point, sq_dist);
481  nearest_facet_recursive(
482  p,
483  nearest_facet, nearest_point, sq_dist,
484  1, 0, mesh_->facets.nb()
485  );
486  return nearest_facet;
487  }
488 
494  index_t nearest_facet(const vec3& p) const {
495  vec3 nearest_point;
496  double sq_dist;
497  return nearest_facet(p, nearest_point, sq_dist);
498  }
499 
520  const vec3& p,
521  index_t& nearest_facet, vec3& nearest_point, double& sq_dist
522  ) const {
523  if(nearest_facet == NO_FACET) {
524  get_nearest_facet_hint(
525  p, nearest_facet, nearest_point, sq_dist
526  );
527  }
528  nearest_facet_recursive(
529  p,
530  nearest_facet, nearest_point, sq_dist,
531  1, 0, mesh_->facets.nb()
532  );
533  }
534 
541  double squared_distance(const vec3& p) const {
542  vec3 nearest_point;
543  double result;
544  nearest_facet(p, nearest_point, result);
545  return result;
546  }
547 
560  const Ray& R,
561  double tmax = Numeric::max_float64(),
562  index_t ignore_f = index_t(-1)
563  ) const;
564 
565 
575  bool ray_nearest_intersection(const Ray& R, Intersection& I) const;
576 
585  bool segment_intersection(const vec3& q1, const vec3& q2) const {
586  return ray_intersection(Ray(q1, q2-q1), 1.0);
587  }
588 
601  const vec3& q1, const vec3& q2, double& t, index_t& f
602  ) const {
603  Ray R(q1, q2-q1);
604  Intersection I;
605  I.t = 1.0;
606  bool result = ray_nearest_intersection(R,I);
607  t = I.t;
608  f = I.f;
609  return result;
610  }
611 
618  const Ray& R,
619  std::function<void(const Intersection&)> action
620  ) const;
621 
622  protected:
623 
638  const vec3& p,
639  index_t& nearest_facet, vec3& nearest_point, double& sq_dist
640  ) const;
641 
660  const vec3& p,
661  index_t& nearest_facet, vec3& nearest_point, double& sq_dist,
662  index_t n, index_t b, index_t e
663  ) const;
664 
681  const Ray& R, const vec3& dirinv, double max_t, index_t ignore_f,
682  index_t n, index_t b, index_t e
683  ) const;
684 
702  const Ray& R, const vec3& dirinv, Intersection& I, index_t ignore_f,
703  index_t n, index_t b, index_t e, index_t coord
704  ) const;
705 
706 
719  const Ray& R, const vec3& dirinv,
720  std::function<void(const Intersection&)> action,
721  index_t n, index_t b, index_t e
722  ) const;
723 
724  };
725 
726  /***********************************************************************/
727 
733  class GEOGRAM_API MeshCellsAABB : public MeshAABB3d {
734  public:
735 
741  static const index_t NO_TET = index_t(-1);
742 
748 
757  MeshCellsAABB(Mesh& M, bool reorder = true);
758 
767  void initialize(Mesh& M, bool reorder = true);
768 
777  index_t containing_tet(const vec3& p) const {
778  geo_debug_assert(mesh_->cells.are_simplices());
779  return containing_tet_recursive(
780  p, 1, 0, mesh_->cells.nb()
781  );
782  }
783 
792  const Box& box_in,
793  std::function<void(index_t)> action
794  ) const {
795  bbox_intersect_recursive(
796  action, box_in, 1, 0, mesh_->cells.nb()
797  );
798  }
799 
808  const vec3& p, std::function<void(index_t)> action
809  ) const {
810  containing_bboxes_recursive(
811  action, p, 1, 0, mesh_->cells.nb()
812  );
813  }
814 
823  std::function<void(index_t, index_t)> action
824  ) const {
825  self_intersect_recursive(
826  action,
827  1, 0, mesh_->cells.nb(),
828  1, 0, mesh_->cells.nb()
829  );
830  }
831 
842  MeshCellsAABB* other,
843  std::function<void(index_t, index_t)> action
844  ) const {
845  other_intersect_recursive(
846  action,
847  1, 0, mesh_->cells.nb(),
848  other,
849  1, 0, other->mesh_->cells.nb()
850  );
851  }
852 
853 
854  protected:
855 
868  const vec3& p,
869  index_t n, index_t b, index_t e
870  ) const;
871 
872 
892  std::function<void(index_t)> action,
893  const vec3& p,
894  index_t node, index_t b, index_t e
895  ) const {
896  geo_debug_assert(e != b);
897 
898  // Prune sub-tree that does not have intersection
899  if(!bboxes_[node].contains(p)) {
900  return;
901  }
902 
903  // Leaf case
904  if(e == b+1) {
905  action(b);
906  return;
907  }
908 
909  // Recursion
910  index_t m = b + (e - b) / 2;
911  index_t node_l = 2 * node;
912  index_t node_r = 2 * node + 1;
913 
914  containing_bboxes_recursive(action, p, node_l, b, m);
915  containing_bboxes_recursive(action, p, node_r, m, e);
916  }
917  };
918 
919 /*******************************************************************/
920 
926  class GEOGRAM_API MeshFacetsAABB2d : public MeshAABB2d {
927  public:
928 
934  static const index_t NO_TRIANGLE = index_t(-1);
935 
941 
950  MeshFacetsAABB2d(Mesh& M, bool reorder = true);
951 
960  void initialize(Mesh& M, bool reorder = true);
961 
970  index_t containing_triangle(const vec2& p) const {
971  geo_debug_assert(mesh_->facets.are_simplices());
972  return containing_triangle_recursive(
973  p, 1, 0, mesh_->facets.nb()
974  );
975  }
976 
985  const Box2d& box_in,
986  std::function<void(index_t)> action
987  ) const {
988  bbox_intersect_recursive(
989  action, box_in, 1, 0, mesh_->facets.nb()
990  );
991  }
992 
1001  const vec2& p, std::function<void(index_t)> action
1002  ) const {
1003  containing_bboxes_recursive(
1004  action, p, 1, 0, mesh_->facets.nb()
1005  );
1006  }
1007 
1016  std::function<void(index_t, index_t)> action
1017  ) const {
1018  self_intersect_recursive(
1019  action,
1020  1, 0, mesh_->facets.nb(),
1021  1, 0, mesh_->facets.nb()
1022  );
1023  }
1024 
1035  MeshFacetsAABB2d* other,
1036  std::function<void(index_t, index_t)> action
1037  ) const {
1038  other_intersect_recursive(
1039  action,
1040  1, 0, mesh_->facets.nb(),
1041  other,
1042  1, 0, other->mesh_->facets.nb()
1043  );
1044  }
1045 
1046 
1047  protected:
1048 
1061  const vec2& p,
1062  index_t n, index_t b, index_t e
1063  ) const;
1064 
1065 
1085  std::function<void(index_t)> action,
1086  const vec2& p,
1087  index_t node, index_t b, index_t e
1088  ) const {
1089  geo_debug_assert(e != b);
1090 
1091  // Prune sub-tree that does not have intersection
1092  if(!bboxes_[node].contains(p)) {
1093  return;
1094  }
1095 
1096  // Leaf case
1097  if(e == b+1) {
1098  action(b);
1099  return;
1100  }
1101 
1102  // Recursion
1103  index_t m = b + (e - b) / 2;
1104  index_t node_l = 2 * node;
1105  index_t node_r = 2 * node + 1;
1106 
1107  containing_bboxes_recursive(action, p, node_l, b, m);
1108  containing_bboxes_recursive(action, p, node_r, m, e);
1109  }
1110  };
1111 
1112 
1113 }
1114 
1115 #endif
#define geo_debug_assert(x)
Verifies that a condition is met.
Definition: assert.h:196
Base class for Axis Aligned Bounding Box trees.
Definition: mesh_AABB.h:60
void bbox_intersect_recursive(std::function< void(index_t)> action, const BOX &box, index_t node, index_t b, index_t e) const
Computes all the elements that have a bbox that intersects a given bbox in a sub-tree of the AABB tre...
Definition: mesh_AABB.h:98
void self_intersect_recursive(std::function< void(index_t, index_t)> action, index_t node1, index_t b1, index_t e1, index_t node2, index_t b2, index_t e2) const
Computes all the pairs of intersecting elements for two sub-trees of the AABB tree.
Definition: mesh_AABB.h:145
void other_intersect_recursive(std::function< void(index_t, index_t)> action, index_t node1, index_t b1, index_t e1, const AABB< BOX > *other, index_t node2, index_t b2, index_t e2) const
Computes all the pairs of intersecting elements for two sub-trees of two AABB trees.
Definition: mesh_AABB.h:214
static index_t max_node_index(index_t node_index, index_t b, index_t e)
Computes the maximum node index in a subtree.
Definition: mesh_AABB.h:269
void initialize(index_t nb, std::function< void(BOX &, index_t)> get_bbox)
Initializes this AABB.
Definition: mesh_AABB.h:70
void init_bboxes_recursive(index_t node_index, index_t b, index_t e, std::function< void(BOX &, index_t)> get_bbox)
Computes the hierarchy of bounding boxes recursively.
Definition: mesh_AABB.h:293
Axis-aligned bounding box.
Definition: geometry.h:732
Axis-aligned bounding box.
Definition: geometry.h:669
Base class for Axis Aligned Bounding Box trees of mesh elements with 2d boxes.
Definition: mesh_AABB.h:330
const Mesh * mesh() const
Gets the mesh.
Definition: mesh_AABB.h:342
MeshAABB2d()
MeshAABB2d constructor.
Definition: mesh_AABB.h:335
Base class for Axis Aligned Bounding Box trees of mesh elements with 3d boxes.
Definition: mesh_AABB.h:356
const Mesh * mesh() const
Gets the mesh.
Definition: mesh_AABB.h:368
MeshAABB3d()
MeshAABB3d constructor.
Definition: mesh_AABB.h:361
Axis Aligned Bounding Box tree of mesh cells.
Definition: mesh_AABB.h:733
void compute_cell_bbox_intersections(std::function< void(index_t, index_t)> action) const
Computes all the pairs of intersecting cells.
Definition: mesh_AABB.h:822
void containing_boxes(const vec3 &p, std::function< void(index_t)> action) const
Finds all the cells such that their bounding box contain a point.
Definition: mesh_AABB.h:807
MeshCellsAABB(Mesh &M, bool reorder=true)
Creates the Axis Aligned Bounding Boxes tree.
MeshCellsAABB()
MeshCellsAABB constructor.
index_t containing_tet_recursive(const vec3 &p, index_t n, index_t b, index_t e) const
The recursive function used by the implementation of containing_tet().
void initialize(Mesh &M, bool reorder=true)
Initializes the Axis Aligned Bounding Boxes tree.
index_t containing_tet(const vec3 &p) const
Finds the index of a tetrahedron that contains a query point.
Definition: mesh_AABB.h:777
void compute_bbox_cell_bbox_intersections(const Box &box_in, std::function< void(index_t)> action) const
Computes all the intersections between a given box and the bounding boxes of all the cells.
Definition: mesh_AABB.h:791
void containing_bboxes_recursive(std::function< void(index_t)> action, const vec3 &p, index_t node, index_t b, index_t e) const
Computes all the cells that have a bbox that contain a given point in a sub-tree of the AABB tree.
Definition: mesh_AABB.h:891
void compute_other_cell_bbox_intersections(MeshCellsAABB *other, std::function< void(index_t, index_t)> action) const
Computes all the pairs of intersecting cells between this AABB and another one.
Definition: mesh_AABB.h:841
Axis Aligned Bounding Box tree of mesh facets in 2D.
Definition: mesh_AABB.h:926
void compute_facet_bbox_intersections(std::function< void(index_t, index_t)> action) const
Computes all the pairs of intersecting facets.
Definition: mesh_AABB.h:1015
void compute_other_cell_bbox_intersections(MeshFacetsAABB2d *other, std::function< void(index_t, index_t)> action) const
Computes all the pairs of intersecting cells between this AABB and another one.
Definition: mesh_AABB.h:1034
index_t containing_triangle_recursive(const vec2 &p, index_t n, index_t b, index_t e) const
The recursive function used by the implementation of containing_triangle().
MeshFacetsAABB2d()
MeshFacetsAABB2d constructor.
void compute_bbox_cell_bbox_intersections(const Box2d &box_in, std::function< void(index_t)> action) const
Computes all the intersections between a given box and the bounding boxes of all the facets.
Definition: mesh_AABB.h:984
void containing_boxes(const vec2 &p, std::function< void(index_t)> action) const
Finds all the cells such that their bounding box contain a point.
Definition: mesh_AABB.h:1000
MeshFacetsAABB2d(Mesh &M, bool reorder=true)
Creates the Axis Aligned Bounding Boxes tree.
void containing_bboxes_recursive(std::function< void(index_t)> action, const vec2 &p, index_t node, index_t b, index_t e) const
Computes all the cells that have a bbox that contain a given point in a sub-tree of the AABB tree.
Definition: mesh_AABB.h:1084
void initialize(Mesh &M, bool reorder=true)
Initializes the Axis Aligned Bounding Boxes tree.
index_t containing_triangle(const vec2 &p) const
Finds the index of a facet that contains a query point.
Definition: mesh_AABB.h:970
Axis Aligned Bounding Box tree of mesh facets in 3D.
Definition: mesh_AABB.h:383
MeshFacetsAABB(Mesh &M, bool reorder=true)
Creates the Axis Aligned Bounding Boxes tree.
bool segment_intersection(const vec3 &q1, const vec3 &q2) const
Tests whether this surface mesh has an intersection with a segment.
Definition: mesh_AABB.h:585
bool ray_intersection_recursive(const Ray &R, const vec3 &dirinv, double max_t, index_t ignore_f, index_t n, index_t b, index_t e) const
The recursive function used by the implementation of ray_intersection()
void nearest_facet_with_hint(const vec3 &p, index_t &nearest_facet, vec3 &nearest_point, double &sq_dist) const
Computes the nearest point and nearest facet from a query point, using user-specified hint.
Definition: mesh_AABB.h:519
void ray_nearest_intersection_recursive(const Ray &R, const vec3 &dirinv, Intersection &I, index_t ignore_f, index_t n, index_t b, index_t e, index_t coord) const
The recursive function used by the implementation of ray_nearest_intersection()
double squared_distance(const vec3 &p) const
Computes the distance between an arbitrary 3d query point and the surface.
Definition: mesh_AABB.h:541
bool ray_nearest_intersection(const Ray &R, Intersection &I) const
Computes the nearest intersection along a ray.
MeshFacetsAABB()
MeshFacetsAABB constructor.
index_t nearest_facet(const vec3 &p, vec3 &nearest_point, double &sq_dist) const
Finds the nearest facet from an arbitrary 3d query point.
Definition: mesh_AABB.h:476
void compute_bbox_facet_bbox_intersections(const Box &box_in, std::function< void(index_t)> action) const
Computes all the intersections between a given box and the bounding boxes of all the facets.
Definition: mesh_AABB.h:461
void get_nearest_facet_hint(const vec3 &p, index_t &nearest_facet, vec3 &nearest_point, double &sq_dist) const
Computes a reasonable initialization for nearest facet search.
void ray_all_intersections_recursive(const Ray &R, const vec3 &dirinv, std::function< void(const Intersection &)> action, index_t n, index_t b, index_t e) const
The function used to implement ray_all_intersections()
void ray_all_intersections(const Ray &R, std::function< void(const Intersection &)> action) const
Calls a user function for all ray-facet intersection.
void compute_facet_bbox_intersections(std::function< void(index_t, index_t)> action) const
Computes all the pairs of intersecting facets.
Definition: mesh_AABB.h:444
void nearest_facet_recursive(const vec3 &p, index_t &nearest_facet, vec3 &nearest_point, double &sq_dist, index_t n, index_t b, index_t e) const
The recursive function used by the implementation of nearest_facet().
bool ray_intersection(const Ray &R, double tmax=Numeric::max_float64(), index_t ignore_f=index_t(-1)) const
Tests whether there exists an intersection between a ray and the mesh.
index_t nearest_facet(const vec3 &p) const
Finds the nearest facet from an arbitrary 3d query point.
Definition: mesh_AABB.h:494
void initialize(Mesh &M, bool reorder=true)
Initializes the Axis Aligned Bounding Boxes tree.
bool segment_nearest_intersection(const vec3 &q1, const vec3 &q2, double &t, index_t &f) const
Finds the intersection between a segment and a surface that is nearest to the first extremity of the ...
Definition: mesh_AABB.h:600
index_t nb() const
Gets the number of (sub-)elements.
Definition: mesh.h:89
Represents a mesh.
Definition: mesh.h:2701
index_t size() const
Gets the number of elements.
Definition: memory.h:674
Common include file, providing basic definitions. Should be included before anything else by all head...
Geometric functions in 2d and 3d.
The class that represents a mesh.
float64 max_float64()
Gets 64 bits float maximum positive value.
Definition: numeric.h:172
Global Vorpaline namespace.
Definition: basic.h:55
bool bboxes_overlap(const Box &B1, const Box &B2)
Tests whether two Boxes have a non-empty intersection.
Definition: geometry.h:701
void get_bbox(const Mesh &M, double *xyzmin, double *xyzmax)
Gets the bounding box of a mesh.
void bbox_union(Box &target, const Box &B1, const Box &B2)
Computes the smallest Box that encloses two Boxes.
Definition: geometry.h:720
geo_index_t index_t
The type for storing and manipulating indices.
Definition: numeric.h:329
Stores all the information related with a ray-facet intersection.
Definition: mesh_AABB.h:390
A Ray, in parametric form.
Definition: geometry.h:971