Geogram Version 1.9.7
A programming library of geometric algorithms
Loading...
Searching...
No Matches
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
50#include <geogram/mesh/mesh.h>
52
53namespace GEO {
54
63 AABB_INPLACE, AABB_INDIRECT, AABB_NOREORDER
64 };
65
66
71 template <class BOX> class AABB {
72
73 protected:
74
82 index_t nb, std::function<void(BOX&, index_t)> get_bbox
83 ) {
84 nb_ = nb;
85 bboxes_.resize(max_node_index(1, 0, nb) + 1);
86 // +1 because size == max_index + 1 !!!
88 }
89
90
110 std::function<void(index_t)> action, const BOX& box,
111 index_t node, index_t b, index_t e
112 ) const {
113 geo_debug_assert(e != b);
114
115 // Prune sub-tree that does not have intersection
116 if(!bboxes_overlap(box, bboxes_[node])) {
117 return;
118 }
119
120 // Leaf case
121 if(e == b+1) {
122 action(element_in_leaf(b));
123 return;
124 }
125
126 // Recursion
127 index_t m = b + (e - b) / 2;
128 index_t node_l = 2 * node;
129 index_t node_r = 2 * node + 1;
130
131 bbox_intersect_recursive(action, box, node_l, b, m);
132 bbox_intersect_recursive(action, box, node_r, m, e);
133 }
134
157 std::function<void(index_t,index_t)> action,
158 index_t node1, index_t b1, index_t e1,
159 index_t node2, index_t b2, index_t e2
160 ) const {
161 geo_debug_assert(e1 != b1);
162 geo_debug_assert(e2 != b2);
163
164 // Since we are intersecting the AABBTree with *itself*,
165 // we can prune half of the cases by skipping the test
166 // whenever node2's facet index interval is greated than
167 // node1's facet index interval.
168 if(e2 <= b1) {
169 return;
170 }
171
172 // The acceleration is here:
173 if(!bboxes_overlap(bboxes_[node1], bboxes_[node2])) {
174 return;
175 }
176
177 // Simple case: leaf - leaf intersection.
178 if(b1 + 1 == e1 && b2 + 1 == e2) {
179 action(element_in_leaf(b1), element_in_leaf(b2));
180 return;
181 }
182
183 // If node2 has more elements than node1, then
184 // intersect node2's two children with node1
185 // else
186 // intersect node1's two children with node2
187 if(e2 - b2 > e1 - b1) {
188 index_t m2 = b2 + (e2 - b2) / 2;
189 index_t node2_l = 2 * node2;
190 index_t node2_r = 2 * node2 + 1;
191 self_intersect_recursive(action, node1, b1, e1, node2_l, b2, m2);
192 self_intersect_recursive(action, node1, b1, e1, node2_r, m2, e2);
193 } else {
194 index_t m1 = b1 + (e1 - b1) / 2;
195 index_t node1_l = 2 * node1;
196 index_t node1_r = 2 * node1 + 1;
197 self_intersect_recursive(action, node1_l, b1, m1, node2, b2, e2);
198 self_intersect_recursive(action, node1_r, m1, e1, node2, b2, e2);
199 }
200 }
201
202
226 std::function<void(index_t,index_t)> action,
227 index_t node1, index_t b1, index_t e1,
228 const AABB<BOX>* other,
229 index_t node2, index_t b2, index_t e2
230 ) const {
231 geo_debug_assert(e1 != b1);
232 geo_debug_assert(e2 != b2);
233
234 // The acceleration is here:
235 if(!bboxes_overlap(bboxes_[node1], other->bboxes_[node2])) {
236 return;
237 }
238
239 // Simple case: leaf - leaf intersection.
240 if(b1 + 1 == e1 && b2 + 1 == e2) {
241 action(element_in_leaf(b1), element_in_leaf(b2));
242 return;
243 }
244
245 // If node2 has more elements than node1, then
246 // intersect node2's two children with node1
247 // else
248 // intersect node1's two children with node2
249 if(e2 - b2 > e1 - b1) {
250 index_t m2 = b2 + (e2 - b2) / 2;
251 index_t node2_l = 2 * node2;
252 index_t node2_r = 2 * node2 + 1;
254 action, node1, b1, e1, other, node2_l, b2, m2
255 );
257 action, node1, b1, e1, other, node2_r, m2, e2
258 );
259 } else {
260 index_t m1 = b1 + (e1 - b1) / 2;
261 index_t node1_l = 2 * node1;
262 index_t node1_r = 2 * node1 + 1;
264 action, node1_l, b1, m1, other, node2, b2, e2
265 );
267 action, node1_r, m1, e1, other, node2, b2, e2
268 );
269 }
270 }
271
272
280 static index_t max_node_index(index_t node_index, index_t b, index_t e) {
281 geo_debug_assert(e > b);
282 if(b + 1 == e) {
283 return node_index;
284 }
285 index_t m = b + (e - b) / 2;
286 index_t childl = 2 * node_index;
287 index_t childr = 2 * node_index + 1;
288 return std::max(
289 max_node_index(childl, b, m),
290 max_node_index(childr, m, e)
291 );
292 }
293
305 index_t node_index,
306 index_t b, index_t e,
307 std::function<void(BOX&, index_t)> get_bbox
308 ) {
309 geo_debug_assert(node_index < bboxes_.size());
310 geo_debug_assert(b != e);
311 if(b + 1 == e) {
312 get_bbox(bboxes_[node_index], element_in_leaf(b));
313 return;
314 }
315 index_t m = b + (e - b) / 2;
316 index_t childl = 2 * node_index;
317 index_t childr = 2 * node_index + 1;
318 geo_debug_assert(childl < bboxes_.size());
319 geo_debug_assert(childr < bboxes_.size());
320 init_bboxes_recursive(childl, b, m, get_bbox);
321 init_bboxes_recursive(childr, m, e, get_bbox);
322 geo_debug_assert(childl < bboxes_.size());
323 geo_debug_assert(childr < bboxes_.size());
324 bbox_union(bboxes_[node_index], bboxes_[childl], bboxes_[childr]);
325 }
326
334 bool indirect() const {
335 return (reorder_.size() != 0);
336 }
337
346 return indirect() ? reorder_[i] : i;
347 }
348
349 protected:
350 index_t nb_;
351 vector<BOX> bboxes_;
353 };
354
355 typedef AABB<Box2d> AABB2d;
356 typedef AABB<Box3d> AABB3d;
357
358 /**************************************************************/
359
364 class GEOGRAM_API MeshAABB2d : public AABB2d {
365 public:
369 MeshAABB2d() : mesh_(nullptr) {
370 }
371
376 const Mesh* mesh() const {
377 return mesh_;
378 }
379
380 protected:
381 Mesh* mesh_;
382 };
383
384 /**************************************************************/
385
390 class GEOGRAM_API MeshAABB3d : public AABB3d {
391 public:
395 MeshAABB3d() : mesh_(nullptr) {
396 }
397
402 const Mesh* mesh() const {
403 return mesh_;
404 }
405
406 protected:
407 Mesh* mesh_;
408 };
409
410 /**************************************************************/
411
417 class GEOGRAM_API MeshFacetsAABB : public MeshAABB3d {
418 public:
419
425 Intersection() :
426 t(Numeric::max_float64()),
427 f(NO_INDEX),
428 i(NO_INDEX), j(NO_INDEX), k(NO_INDEX)
429 {
430 }
432 double t;
435 index_t i,j,k;
436 double u,v;
437 };
438
439
445 }
446
461 initialize(M, reorder_mode);
462 }
463
470 initialize(const_cast<Mesh&>(M), AABB_INDIRECT);
471 }
472
486 void initialize(Mesh& M, AABBReorderMode reorder_mode = AABB_INDIRECT);
487
488#ifndef GOMGEN
489 [[deprecated("use MeshFacetsAABB(Mesh&,AABBReorderMode) instead")]]
490 MeshFacetsAABB(Mesh& M, bool reorder) {
491 initialize(M, reorder ? AABB_INPLACE : AABB_NOREORDER);
492 }
493
494 [[deprecated("use initialize(Mesh&,AABBReorderMode) instead")]]
495 void initialize(Mesh& M, bool reorder) {
496 initialize(M, reorder ? AABB_INPLACE : AABB_NOREORDER);
497 }
498#endif
507 std::function<void(index_t, index_t)> action
508 ) const {
509 self_intersect_recursive(
510 action,
511 1, 0, mesh_->facets.nb(),
512 1, 0, mesh_->facets.nb()
513 );
514 }
515
524 const Box& box_in, std::function<void(index_t)> action
525 ) const {
526 bbox_intersect_recursive(
527 action, box_in, 1, 0, mesh_->facets.nb()
528 );
529 }
530
539 const vec3& p, vec3& nearest_point, double& sq_dist
540 ) const {
541 index_t nearest_facet;
542 get_nearest_facet_hint(p, nearest_facet, nearest_point, sq_dist);
543 nearest_facet_recursive(
544 p,
545 nearest_facet, nearest_point, sq_dist,
546 1, 0, mesh_->facets.nb()
547 );
548 return nearest_facet;
549 }
550
556 index_t nearest_facet(const vec3& p) const {
557 vec3 nearest_point;
558 double sq_dist;
559 return nearest_facet(p, nearest_point, sq_dist);
560 }
561
562
575 const vec3& p, vec3& nearest_point, double& sq_dist,
576 std::function<bool(index_t)> filter
577 ) const {
578 index_t nearest_facet = NO_INDEX;
579 sq_dist = Numeric::max_float64();
580 nearest_facet_recursive_filtered(
581 p,
582 nearest_facet, nearest_point, sq_dist,
583 1, 0, mesh_->facets.nb(),
584 filter
585 );
586 return nearest_facet;
587 }
588
609 const vec3& p,
610 index_t& nearest_facet, vec3& nearest_point, double& sq_dist
611 ) const {
612 if(nearest_facet == NO_FACET) {
613 get_nearest_facet_hint(
614 p, nearest_facet, nearest_point, sq_dist
615 );
616 }
617 nearest_facet_recursive(
618 p,
619 nearest_facet, nearest_point, sq_dist,
620 1, 0, mesh_->facets.nb()
621 );
622 }
623
630 double squared_distance(const vec3& p) const {
631 vec3 nearest_point;
632 double result;
633 nearest_facet(p, nearest_point, result);
634 return result;
635 }
636
649 const Ray& R,
650 double tmax = Numeric::max_float64(),
651 index_t ignore_f = NO_INDEX
652 ) const;
653
654
664 bool ray_nearest_intersection(const Ray& R, Intersection& I) const;
665
674 bool segment_intersection(const vec3& q1, const vec3& q2) const {
675 return ray_intersection(Ray(q1, q2-q1), 1.0);
676 }
677
690 const vec3& q1, const vec3& q2, double& t, index_t& f
691 ) const {
692 Ray R(q1, q2-q1);
693 Intersection I;
694 I.t = 1.0;
695 bool result = ray_nearest_intersection(R,I);
696 t = I.t;
697 f = I.f;
698 return result;
699 }
700
707 const Ray& R,
708 std::function<void(const Intersection&)> action
709 ) const;
710
711
720 bool contains(const vec3& p) const;
721
722 protected:
723
738 const vec3& p,
739 index_t& nearest_facet, vec3& nearest_point, double& sq_dist
740 ) const;
741
760 const vec3& p,
761 index_t& nearest_facet, vec3& nearest_point, double& sq_dist,
762 index_t n, index_t b, index_t e
763 ) const;
764
765
783 const vec3& p,
784 index_t& nearest_facet, vec3& nearest_point, double& sq_dist,
785 index_t n, index_t b, index_t e,
786 std::function<bool(index_t)> filter
787 ) const;
788
805 const Ray& R, const vec3& dirinv, double max_t, index_t ignore_f,
806 index_t n, index_t b, index_t e
807 ) const;
808
826 const Ray& R, const vec3& dirinv, Intersection& I, index_t ignore_f,
827 index_t n, index_t b, index_t e, index_t coord
828 ) const;
829
830
843 const Ray& R, const vec3& dirinv,
844 std::function<void(const Intersection&)> action,
845 index_t n, index_t b, index_t e
846 ) const;
847
848 };
849
850 /***********************************************************************/
851
857 class GEOGRAM_API MeshCellsAABB : public MeshAABB3d {
858 public:
859
865 static constexpr index_t NO_TET = NO_INDEX;
866
872 }
873
886 initialize(M, reorder_mode);
887 }
888
894 MeshCellsAABB(const Mesh& M) {
895 initialize(const_cast<Mesh&>(M), AABB_INDIRECT);
896 }
897
909 void initialize(Mesh& M, AABBReorderMode reorder_mode = AABB_INDIRECT);
910
911#ifndef GOMGEN
912 [[deprecated("use MeshCellsAABB(Mesh&,AABBReorderMode) instead")]]
913 MeshCellsAABB(Mesh& M, bool reorder) {
914 initialize(M, reorder ? AABB_INPLACE : AABB_NOREORDER);
915 }
916
917 [[deprecated("use initialize(Mesh&,AABBReorderMode) instead")]]
918 void initialize(Mesh& M, bool reorder) {
919 initialize(M, reorder ? AABB_INPLACE : AABB_NOREORDER);
920 }
921#endif
922
931 index_t containing_tet(const vec3& p) const {
932 geo_debug_assert(mesh_->cells.are_simplices());
933 return containing_tet_recursive(
934 p, 1, 0, mesh_->cells.nb()
935 );
936 }
937
946 const Box& box_in,
947 std::function<void(index_t)> action
948 ) const {
949 bbox_intersect_recursive(
950 action, box_in, 1, 0, mesh_->cells.nb()
951 );
952 }
953
962 const vec3& p, std::function<void(index_t)> action
963 ) const {
964 containing_bboxes_recursive(
965 action, p, 1, 0, mesh_->cells.nb()
966 );
967 }
968
977 std::function<void(index_t, index_t)> action
978 ) const {
979 self_intersect_recursive(
980 action,
981 1, 0, mesh_->cells.nb(),
982 1, 0, mesh_->cells.nb()
983 );
984 }
985
996 MeshCellsAABB* other,
997 std::function<void(index_t, index_t)> action
998 ) const {
999 other_intersect_recursive(
1000 action,
1001 1, 0, mesh_->cells.nb(),
1002 other,
1003 1, 0, other->mesh_->cells.nb()
1004 );
1005 }
1006
1007
1008 protected:
1009
1022 const vec3& p,
1023 index_t n, index_t b, index_t e
1024 ) const;
1025
1026
1046 std::function<void(index_t)> action,
1047 const vec3& p,
1048 index_t node, index_t b, index_t e
1049 ) const {
1050 geo_debug_assert(e != b);
1051
1052 // Prune sub-tree that does not have intersection
1053 if(!bboxes_[node].contains(p)) {
1054 return;
1055 }
1056
1057 // Leaf case
1058 if(e == b+1) {
1059 action(element_in_leaf(b));
1060 return;
1061 }
1062
1063 // Recursion
1064 index_t m = b + (e - b) / 2;
1065 index_t node_l = 2 * node;
1066 index_t node_r = 2 * node + 1;
1067
1068 containing_bboxes_recursive(action, p, node_l, b, m);
1069 containing_bboxes_recursive(action, p, node_r, m, e);
1070 }
1071 };
1072
1073/*******************************************************************/
1074
1080 class GEOGRAM_API MeshFacetsAABB2d : public MeshAABB2d {
1081 public:
1082
1088 static constexpr index_t NO_TRIANGLE = NO_INDEX;
1089
1095
1104 MeshFacetsAABB2d(Mesh& M, bool reorder = true);
1105
1114 void initialize(Mesh& M, bool reorder = true);
1115
1125 geo_debug_assert(mesh_->facets.are_simplices());
1126 return containing_triangle_recursive(
1127 p, 1, 0, mesh_->facets.nb()
1128 );
1129 }
1130
1139 const Box2d& box_in,
1140 std::function<void(index_t)> action
1141 ) const {
1142 bbox_intersect_recursive(
1143 action, box_in, 1, 0, mesh_->facets.nb()
1144 );
1145 }
1146
1155 const vec2& p, std::function<void(index_t)> action
1156 ) const {
1157 containing_bboxes_recursive(
1158 action, p, 1, 0, mesh_->facets.nb()
1159 );
1160 }
1161
1170 std::function<void(index_t, index_t)> action
1171 ) const {
1172 self_intersect_recursive(
1173 action,
1174 1, 0, mesh_->facets.nb(),
1175 1, 0, mesh_->facets.nb()
1176 );
1177 }
1178
1189 MeshFacetsAABB2d* other,
1190 std::function<void(index_t, index_t)> action
1191 ) const {
1192 other_intersect_recursive(
1193 action,
1194 1, 0, mesh_->facets.nb(),
1195 other,
1196 1, 0, other->mesh_->facets.nb()
1197 );
1198 }
1199
1200
1201 protected:
1202
1215 const vec2& p,
1216 index_t n, index_t b, index_t e
1217 ) const;
1218
1219
1239 std::function<void(index_t)> action,
1240 const vec2& p,
1241 index_t node, index_t b, index_t e
1242 ) const {
1243 geo_debug_assert(e != b);
1244
1245 // Prune sub-tree that does not have intersection
1246 if(!bboxes_[node].contains(p)) {
1247 return;
1248 }
1249
1250 // Leaf case
1251 if(e == b+1) {
1252 action(element_in_leaf(b));
1253 return;
1254 }
1255
1256 // Recursion
1257 index_t m = b + (e - b) / 2;
1258 index_t node_l = 2 * node;
1259 index_t node_r = 2 * node + 1;
1260
1261 containing_bboxes_recursive(action, p, node_l, b, m);
1262 containing_bboxes_recursive(action, p, node_r, m, e);
1263 }
1264 };
1265
1266
1267}
1268
1269#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:71
bool indirect() const
Tests whether this AABB is indirect or in-place.
Definition mesh_AABB.h:334
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:109
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:156
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:225
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:280
vector< index_t > reorder_
Definition mesh_AABB.h:352
void initialize(index_t nb, std::function< void(BOX &, index_t)> get_bbox)
Initializes this AABB.
Definition mesh_AABB.h:81
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:304
index_t element_in_leaf(index_t i) const
Gets the element stored in a leaf node.
Definition mesh_AABB.h:345
Axis-aligned bounding box.
Definition geometry.h:752
Axis-aligned bounding box.
Definition geometry.h:689
Base class for Axis Aligned Bounding Box trees of mesh elements with 2d boxes.
Definition mesh_AABB.h:364
const Mesh * mesh() const
Gets the mesh.
Definition mesh_AABB.h:376
MeshAABB2d()
MeshAABB2d constructor.
Definition mesh_AABB.h:369
Base class for Axis Aligned Bounding Box trees of mesh elements with 3d boxes.
Definition mesh_AABB.h:390
const Mesh * mesh() const
Gets the mesh.
Definition mesh_AABB.h:402
MeshAABB3d()
MeshAABB3d constructor.
Definition mesh_AABB.h:395
Axis Aligned Bounding Box tree of mesh cells.
Definition mesh_AABB.h:857
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:976
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:961
MeshCellsAABB()
MeshCellsAABB constructor.
Definition mesh_AABB.h:871
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().
index_t containing_tet(const vec3 &p) const
Finds the index of a tetrahedron that contains a query point.
Definition mesh_AABB.h:931
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:945
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:1045
MeshCellsAABB(Mesh &M, AABBReorderMode reorder_mode)
Creates an Axis Aligned Bounding Boxes tree for mesh cells.
Definition mesh_AABB.h:885
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:995
MeshCellsAABB(const Mesh &M)
Creates an Axis Aligned Bounding Boxes tree for mesh cells.
Definition mesh_AABB.h:894
void initialize(Mesh &M, AABBReorderMode reorder_mode=AABB_INDIRECT)
Initializes the Axis Aligned Bounding Boxes tree.
Axis Aligned Bounding Box tree of mesh facets in 2D.
Definition mesh_AABB.h:1080
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:1169
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:1188
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:1138
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:1154
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:1238
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:1124
Axis Aligned Bounding Box tree of mesh facets in 3D.
Definition mesh_AABB.h:417
bool ray_intersection(const Ray &R, double tmax=Numeric::max_float64(), index_t ignore_f=NO_INDEX) const
Tests whether there exists an intersection between a ray and the mesh.
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:674
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 initialize(Mesh &M, AABBReorderMode reorder_mode=AABB_INDIRECT)
Initializes the Axis Aligned Bounding Boxes tree.
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:608
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()
void nearest_facet_recursive_filtered(const vec3 &p, index_t &nearest_facet, vec3 &nearest_point, double &sq_dist, index_t n, index_t b, index_t e, std::function< bool(index_t)> filter) const
The recursive function used by the implementation of nearest_facet_filtered().
double squared_distance(const vec3 &p) const
Computes the distance between an arbitrary 3d query point and the surface.
Definition mesh_AABB.h:630
bool ray_nearest_intersection(const Ray &R, Intersection &I) const
Computes the nearest intersection along a ray.
MeshFacetsAABB(Mesh &M, AABBReorderMode reorder_mode)
Creates an Axis Aligned Bounding Boxes tree for facets.
Definition mesh_AABB.h:460
MeshFacetsAABB()
MeshFacetsAABB constructor.
Definition mesh_AABB.h:444
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:538
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:523
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:506
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 contains(const vec3 &p) const
Tests whether a closed surface contains a point.
MeshFacetsAABB(const Mesh &M)
Creates an Axis Aligned Bounding Boxes tree for facets.
Definition mesh_AABB.h:469
index_t nearest_facet_filtered(const vec3 &p, vec3 &nearest_point, double &sq_dist, std::function< bool(index_t)> filter) const
Finds the nearest facet from an arbitrary 3d query point taken into account only certain facets.
Definition mesh_AABB.h:574
index_t nearest_facet(const vec3 &p) const
Finds the nearest facet from an arbitrary 3d query point.
Definition mesh_AABB.h:556
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:689
index_t nb() const
Gets the number of (sub-)elements.
Definition mesh.h:98
Represents a mesh.
Definition mesh.h:3132
Vector with aligned memory allocation.
Definition memory.h:660
index_t size() const
Gets the number of elements.
Definition memory.h:706
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.
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:721
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:740
AABBReorderMode
Axis Aligned Bounding Box Tree ordering mode.
Definition mesh_AABB.h:62
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:424
A Ray, in parametric form.
Definition geometry.h:991