Geogram  Version 1.9.1-rc
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(index_t nb, std::function<void(BOX&, index_t)> get_bbox) {
71  nb_ = nb;
72  bboxes_.resize(max_node_index(1, 0, nb) + 1);
73  // +1 because size == max_index + 1 !!!
74  init_bboxes_recursive(1, 0, nb_, get_bbox);
75  }
76 
77 
97  std::function<void(index_t)> action, const BOX& box,
98  index_t node, index_t b, index_t e
99  ) const {
100  geo_debug_assert(e != b);
101 
102  // Prune sub-tree that does not have intersection
103  if(!bboxes_overlap(box, bboxes_[node])) {
104  return;
105  }
106 
107  // Leaf case
108  if(e == b+1) {
109  action(b);
110  return;
111  }
112 
113  // Recursion
114  index_t m = b + (e - b) / 2;
115  index_t node_l = 2 * node;
116  index_t node_r = 2 * node + 1;
117 
118  bbox_intersect_recursive(action, box, node_l, b, m);
119  bbox_intersect_recursive(action, box, node_r, m, e);
120  }
121 
144  std::function<void(index_t,index_t)> action,
145  index_t node1, index_t b1, index_t e1,
146  index_t node2, index_t b2, index_t e2
147  ) const {
148  geo_debug_assert(e1 != b1);
149  geo_debug_assert(e2 != b2);
150 
151  // Since we are intersecting the AABBTree with *itself*,
152  // we can prune half of the cases by skipping the test
153  // whenever node2's facet index interval is greated than
154  // node1's facet index interval.
155  if(e2 <= b1) {
156  return;
157  }
158 
159  // The acceleration is here:
160  if(!bboxes_overlap(bboxes_[node1], bboxes_[node2])) {
161  return;
162  }
163 
164  // Simple case: leaf - leaf intersection.
165  if(b1 + 1 == e1 && b2 + 1 == e2) {
166  action(b1, b2);
167  return;
168  }
169 
170  // If node2 has more elements than node1, then
171  // intersect node2's two children with node1
172  // else
173  // intersect node1's two children with node2
174  if(e2 - b2 > e1 - b1) {
175  index_t m2 = b2 + (e2 - b2) / 2;
176  index_t node2_l = 2 * node2;
177  index_t node2_r = 2 * node2 + 1;
178  self_intersect_recursive(action, node1, b1, e1, node2_l, b2, m2);
179  self_intersect_recursive(action, node1, b1, e1, node2_r, m2, e2);
180  } else {
181  index_t m1 = b1 + (e1 - b1) / 2;
182  index_t node1_l = 2 * node1;
183  index_t node1_r = 2 * node1 + 1;
184  self_intersect_recursive(action, node1_l, b1, m1, node2, b2, e2);
185  self_intersect_recursive(action, node1_r, m1, e1, node2, b2, e2);
186  }
187  }
188 
189 
213  std::function<void(index_t,index_t)> action,
214  index_t node1, index_t b1, index_t e1,
215  const AABB<BOX>* other,
216  index_t node2, index_t b2, index_t e2
217  ) const {
218  geo_debug_assert(e1 != b1);
219  geo_debug_assert(e2 != b2);
220 
221  // The acceleration is here:
222  if(!bboxes_overlap(bboxes_[node1], other->bboxes_[node2])) {
223  return;
224  }
225 
226  // Simple case: leaf - leaf intersection.
227  if(b1 + 1 == e1 && b2 + 1 == e2) {
228  action(b1, b2);
229  return;
230  }
231 
232  // If node2 has more elements than node1, then
233  // intersect node2's two children with node1
234  // else
235  // intersect node1's two children with node2
236  if(e2 - b2 > e1 - b1) {
237  index_t m2 = b2 + (e2 - b2) / 2;
238  index_t node2_l = 2 * node2;
239  index_t node2_r = 2 * node2 + 1;
241  action, node1, b1, e1, other, node2_l, b2, m2
242  );
244  action, node1, b1, e1, other, node2_r, m2, e2
245  );
246  } else {
247  index_t m1 = b1 + (e1 - b1) / 2;
248  index_t node1_l = 2 * node1;
249  index_t node1_r = 2 * node1 + 1;
251  action, node1_l, b1, m1, other, node2, b2, e2
252  );
254  action, node1_r, m1, e1, other, node2, b2, e2
255  );
256  }
257  }
258 
259 
267  static index_t max_node_index(index_t node_index, index_t b, index_t e) {
268  geo_debug_assert(e > b);
269  if(b + 1 == e) {
270  return node_index;
271  }
272  index_t m = b + (e - b) / 2;
273  index_t childl = 2 * node_index;
274  index_t childr = 2 * node_index + 1;
275  return std::max(
276  max_node_index(childl, b, m),
277  max_node_index(childr, m, e)
278  );
279  }
280 
292  index_t node_index,
293  index_t b, index_t e,
294  std::function<void(BOX&, index_t)> get_bbox
295  ) {
296  geo_debug_assert(node_index < bboxes_.size());
297  geo_debug_assert(b != e);
298  if(b + 1 == e) {
299  get_bbox(bboxes_[node_index], b);
300  return;
301  }
302  index_t m = b + (e - b) / 2;
303  index_t childl = 2 * node_index;
304  index_t childr = 2 * node_index + 1;
305  geo_debug_assert(childl < bboxes_.size());
306  geo_debug_assert(childr < bboxes_.size());
307  init_bboxes_recursive(childl, b, m, get_bbox);
308  init_bboxes_recursive(childr, m, e, get_bbox);
309  geo_debug_assert(childl < bboxes_.size());
310  geo_debug_assert(childr < bboxes_.size());
311  bbox_union(bboxes_[node_index], bboxes_[childl], bboxes_[childr]);
312  }
313 
314  protected:
315  index_t nb_;
316  vector<BOX> bboxes_;
317  };
318 
319  typedef AABB<Box2d> AABB2d;
320  typedef AABB<Box3d> AABB3d;
321 
322  /**************************************************************/
323 
328  class GEOGRAM_API MeshAABB2d : public AABB2d {
329  public:
333  MeshAABB2d() : mesh_(nullptr) {
334  }
335 
340  const Mesh* mesh() const {
341  return mesh_;
342  }
343 
344  protected:
345  Mesh* mesh_;
346  };
347 
348  /**************************************************************/
349 
354  class GEOGRAM_API MeshAABB3d : public AABB3d {
355  public:
359  MeshAABB3d() : mesh_(nullptr) {
360  }
361 
366  const Mesh* mesh() const {
367  return mesh_;
368  }
369 
370  protected:
371  Mesh* mesh_;
372  };
373 
374  /**************************************************************/
375 
381  class GEOGRAM_API MeshFacetsAABB : public MeshAABB3d {
382  public:
383 
388  struct Intersection {
389  Intersection() :
391  f(index_t(-1)),
392  i(index_t(-1)), j(index_t(-1)), k(index_t(-1))
393  {
394  }
395  vec3 p;
396  double t;
398  vec3 N;
399  index_t i,j,k;
400  double u,v;
401  };
402 
403 
409 
420  void initialize(Mesh& M, bool reorder = true);
421 
422 
433  MeshFacetsAABB(Mesh& M, bool reorder = true);
434 
443  std::function<void(index_t, index_t)> action
444  ) const {
445  self_intersect_recursive(
446  action,
447  1, 0, mesh_->facets.nb(),
448  1, 0, mesh_->facets.nb()
449  );
450  }
451 
460  const Box& box_in, std::function<void(index_t)> action
461  ) const {
462  bbox_intersect_recursive(
463  action, box_in, 1, 0, mesh_->facets.nb()
464  );
465  }
466 
475  const vec3& p, vec3& nearest_point, double& sq_dist
476  ) const {
477  index_t nearest_facet;
478  get_nearest_facet_hint(p, nearest_facet, nearest_point, sq_dist);
479  nearest_facet_recursive(
480  p,
481  nearest_facet, nearest_point, sq_dist,
482  1, 0, mesh_->facets.nb()
483  );
484  return nearest_facet;
485  }
486 
492  index_t nearest_facet(const vec3& p) const {
493  vec3 nearest_point;
494  double sq_dist;
495  return nearest_facet(p, nearest_point, sq_dist);
496  }
497 
518  const vec3& p,
519  index_t& nearest_facet, vec3& nearest_point, double& sq_dist
520  ) const {
521  if(nearest_facet == NO_FACET) {
522  get_nearest_facet_hint(
523  p, nearest_facet, nearest_point, sq_dist
524  );
525  }
526  nearest_facet_recursive(
527  p,
528  nearest_facet, nearest_point, sq_dist,
529  1, 0, mesh_->facets.nb()
530  );
531  }
532 
539  double squared_distance(const vec3& p) const {
540  vec3 nearest_point;
541  double result;
542  nearest_facet(p, nearest_point, result);
543  return result;
544  }
545 
558  const Ray& R,
559  double tmax = Numeric::max_float64(),
560  index_t ignore_f = index_t(-1)
561  ) const;
562 
563 
573  bool ray_nearest_intersection(const Ray& R, Intersection& I) const;
574 
583  bool segment_intersection(const vec3& q1, const vec3& q2) const {
584  return ray_intersection(Ray(q1, q2-q1), 1.0);
585  }
586 
599  const vec3& q1, const vec3& q2, double& t, index_t& f
600  ) const {
601  Ray R(q1, q2-q1);
602  Intersection I;
603  I.t = 1.0;
604  bool result = ray_nearest_intersection(R,I);
605  t = I.t;
606  f = I.f;
607  return result;
608  }
609 
616  const Ray& R,
617  std::function<void(const Intersection&)> action
618  ) const;
619 
620  protected:
621 
636  const vec3& p,
637  index_t& nearest_facet, vec3& nearest_point, double& sq_dist
638  ) const;
639 
658  const vec3& p,
659  index_t& nearest_facet, vec3& nearest_point, double& sq_dist,
660  index_t n, index_t b, index_t e
661  ) const;
662 
679  const Ray& R, const vec3& dirinv, double max_t, index_t ignore_f,
680  index_t n, index_t b, index_t e
681  ) const;
682 
700  const Ray& R, const vec3& dirinv, Intersection& I, index_t ignore_f,
701  index_t n, index_t b, index_t e, index_t coord
702  ) const;
703 
704 
717  const Ray& R, const vec3& dirinv,
718  std::function<void(const Intersection&)> action,
719  index_t n, index_t b, index_t e
720  ) const;
721 
722  };
723 
724  /***********************************************************************/
725 
731  class GEOGRAM_API MeshCellsAABB : public MeshAABB3d {
732  public:
733 
739  static const index_t NO_TET = index_t(-1);
740 
746 
755  MeshCellsAABB(Mesh& M, bool reorder = true);
756 
765  void initialize(Mesh& M, bool reorder = true);
766 
775  index_t containing_tet(const vec3& p) const {
776  geo_debug_assert(mesh_->cells.are_simplices());
777  return containing_tet_recursive(
778  p, 1, 0, mesh_->cells.nb()
779  );
780  }
781 
790  const Box& box_in,
791  std::function<void(index_t)> action
792  ) const {
793  bbox_intersect_recursive(
794  action, box_in, 1, 0, mesh_->cells.nb()
795  );
796  }
797 
806  const vec3& p, std::function<void(index_t)> action
807  ) const {
808  containing_bboxes_recursive(
809  action, p, 1, 0, mesh_->cells.nb()
810  );
811  }
812 
821  std::function<void(index_t, index_t)> action
822  ) const {
823  self_intersect_recursive(
824  action,
825  1, 0, mesh_->cells.nb(),
826  1, 0, mesh_->cells.nb()
827  );
828  }
829 
840  MeshCellsAABB* other,
841  std::function<void(index_t, index_t)> action
842  ) const {
843  other_intersect_recursive(
844  action,
845  1, 0, mesh_->cells.nb(),
846  other,
847  1, 0, other->mesh_->cells.nb()
848  );
849  }
850 
851 
852  protected:
853 
866  const vec3& p,
867  index_t n, index_t b, index_t e
868  ) const;
869 
870 
890  std::function<void(index_t)> action,
891  const vec3& p,
892  index_t node, index_t b, index_t e
893  ) const {
894  geo_debug_assert(e != b);
895 
896  // Prune sub-tree that does not have intersection
897  if(!bboxes_[node].contains(p)) {
898  return;
899  }
900 
901  // Leaf case
902  if(e == b+1) {
903  action(b);
904  return;
905  }
906 
907  // Recursion
908  index_t m = b + (e - b) / 2;
909  index_t node_l = 2 * node;
910  index_t node_r = 2 * node + 1;
911 
912  containing_bboxes_recursive(action, p, node_l, b, m);
913  containing_bboxes_recursive(action, p, node_r, m, e);
914  }
915  };
916 
917 /*******************************************************************/
918 
924  class GEOGRAM_API MeshFacetsAABB2d : public MeshAABB2d {
925  public:
926 
932  static const index_t NO_TRIANGLE = index_t(-1);
933 
939 
948  MeshFacetsAABB2d(Mesh& M, bool reorder = true);
949 
958  void initialize(Mesh& M, bool reorder = true);
959 
968  index_t containing_triangle(const vec2& p) const {
969  geo_debug_assert(mesh_->facets.are_simplices());
970  return containing_triangle_recursive(
971  p, 1, 0, mesh_->facets.nb()
972  );
973  }
974 
983  const Box2d& box_in,
984  std::function<void(index_t)> action
985  ) const {
986  bbox_intersect_recursive(
987  action, box_in, 1, 0, mesh_->facets.nb()
988  );
989  }
990 
999  const vec2& p, std::function<void(index_t)> action
1000  ) const {
1001  containing_bboxes_recursive(
1002  action, p, 1, 0, mesh_->facets.nb()
1003  );
1004  }
1005 
1014  std::function<void(index_t, index_t)> action
1015  ) const {
1016  self_intersect_recursive(
1017  action,
1018  1, 0, mesh_->facets.nb(),
1019  1, 0, mesh_->facets.nb()
1020  );
1021  }
1022 
1033  MeshFacetsAABB2d* other,
1034  std::function<void(index_t, index_t)> action
1035  ) const {
1036  other_intersect_recursive(
1037  action,
1038  1, 0, mesh_->facets.nb(),
1039  other,
1040  1, 0, other->mesh_->facets.nb()
1041  );
1042  }
1043 
1044 
1045  protected:
1046 
1059  const vec2& p,
1060  index_t n, index_t b, index_t e
1061  ) const;
1062 
1063 
1083  std::function<void(index_t)> action,
1084  const vec2& p,
1085  index_t node, index_t b, index_t e
1086  ) const {
1087  geo_debug_assert(e != b);
1088 
1089  // Prune sub-tree that does not have intersection
1090  if(!bboxes_[node].contains(p)) {
1091  return;
1092  }
1093 
1094  // Leaf case
1095  if(e == b+1) {
1096  action(b);
1097  return;
1098  }
1099 
1100  // Recursion
1101  index_t m = b + (e - b) / 2;
1102  index_t node_l = 2 * node;
1103  index_t node_r = 2 * node + 1;
1104 
1105  containing_bboxes_recursive(action, p, node_l, b, m);
1106  containing_bboxes_recursive(action, p, node_r, m, e);
1107  }
1108  };
1109 
1110 
1111 }
1112 
1113 #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:96
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:143
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:212
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:267
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:291
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:328
const Mesh * mesh() const
Gets the mesh.
Definition: mesh_AABB.h:340
MeshAABB2d()
MeshAABB2d constructor.
Definition: mesh_AABB.h:333
Base class for Axis Aligned Bounding Box trees of mesh elements with 3d boxes.
Definition: mesh_AABB.h:354
const Mesh * mesh() const
Gets the mesh.
Definition: mesh_AABB.h:366
MeshAABB3d()
MeshAABB3d constructor.
Definition: mesh_AABB.h:359
Axis Aligned Bounding Box tree of mesh cells.
Definition: mesh_AABB.h:731
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:820
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:805
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:775
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:789
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:889
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:839
Axis Aligned Bounding Box tree of mesh facets in 2D.
Definition: mesh_AABB.h:924
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:1013
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:1032
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:982
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:998
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:1082
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:968
Axis Aligned Bounding Box tree of mesh facets in 3D.
Definition: mesh_AABB.h:381
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:583
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:517
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:539
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:474
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:459
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:442
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:492
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:598
index_t nb() const
Gets the number of (sub-)elements.
Definition: mesh.h:89
Represents a mesh.
Definition: mesh.h:2693
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:388
A Ray, in parametric form.
Definition: geometry.h:971