1// Copyright 2009-2021 Intel Corporation
2// SPDX-License-Identifier: Apache-2.0
3
4#include "bvh_collider.h"
5#include "../geometry/triangle_triangle_intersector.h"
6
7namespace embree
8{
9 namespace isa
10 {
11#define CSTAT(x)
12
13 size_t parallel_depth_threshold = 3;
14 CSTAT(std::atomic<size_t> bvh_collide_traversal_steps(0));
15 CSTAT(std::atomic<size_t> bvh_collide_leaf_pairs(0));
16 CSTAT(std::atomic<size_t> bvh_collide_leaf_iterations(0));
17 CSTAT(std::atomic<size_t> bvh_collide_prim_intersections1(0));
18 CSTAT(std::atomic<size_t> bvh_collide_prim_intersections2(0));
19 CSTAT(std::atomic<size_t> bvh_collide_prim_intersections3(0));
20 CSTAT(std::atomic<size_t> bvh_collide_prim_intersections4(0));
21 CSTAT(std::atomic<size_t> bvh_collide_prim_intersections5(0));
22 CSTAT(std::atomic<size_t> bvh_collide_prim_intersections(0));
23
24 struct Collision
25 {
26 __forceinline Collision() {}
27
28 __forceinline Collision (unsigned geomID0, unsigned primID0, unsigned geomID1, unsigned primID1)
29 : geomID0(geomID0), primID0(primID0), geomID1(geomID1), primID1(primID1) {}
30
31 unsigned geomID0;
32 unsigned primID0;
33 unsigned geomID1;
34 unsigned primID1;
35 };
36
37 template<int N>
38 __forceinline size_t overlap(const BBox3fa& box0, const typename BVHN<N>::AABBNode& node1)
39 {
40 const vfloat<N> lower_x = max(vfloat<N>(box0.lower.x),node1.lower_x);
41 const vfloat<N> lower_y = max(vfloat<N>(box0.lower.y),node1.lower_y);
42 const vfloat<N> lower_z = max(vfloat<N>(box0.lower.z),node1.lower_z);
43 const vfloat<N> upper_x = min(vfloat<N>(box0.upper.x),node1.upper_x);
44 const vfloat<N> upper_y = min(vfloat<N>(box0.upper.y),node1.upper_y);
45 const vfloat<N> upper_z = min(vfloat<N>(box0.upper.z),node1.upper_z);
46 return movemask((lower_x <= upper_x) & (lower_y <= upper_y) & (lower_z <= upper_z));
47 }
48
49 template<int N>
50 __forceinline size_t overlap(const BBox3fa& box0, const BBox<Vec3<vfloat<N>>>& box1)
51 {
52 const vfloat<N> lower_x = max(vfloat<N>(box0.lower.x),box1.lower.x);
53 const vfloat<N> lower_y = max(vfloat<N>(box0.lower.y),box1.lower.y);
54 const vfloat<N> lower_z = max(vfloat<N>(box0.lower.z),box1.lower.z);
55 const vfloat<N> upper_x = min(vfloat<N>(box0.upper.x),box1.upper.x);
56 const vfloat<N> upper_y = min(vfloat<N>(box0.upper.y),box1.upper.y);
57 const vfloat<N> upper_z = min(vfloat<N>(box0.upper.z),box1.upper.z);
58 return movemask((lower_x <= upper_x) & (lower_y <= upper_y) & (lower_z <= upper_z));
59 }
60
61 template<int N>
62 __forceinline size_t overlap(const BBox<Vec3<vfloat<N>>>& box0, size_t i, const BBox<Vec3<vfloat<N>>>& box1)
63 {
64 const vfloat<N> lower_x = max(vfloat<N>(box0.lower.x[i]),box1.lower.x);
65 const vfloat<N> lower_y = max(vfloat<N>(box0.lower.y[i]),box1.lower.y);
66 const vfloat<N> lower_z = max(vfloat<N>(box0.lower.z[i]),box1.lower.z);
67 const vfloat<N> upper_x = min(vfloat<N>(box0.upper.x[i]),box1.upper.x);
68 const vfloat<N> upper_y = min(vfloat<N>(box0.upper.y[i]),box1.upper.y);
69 const vfloat<N> upper_z = min(vfloat<N>(box0.upper.z[i]),box1.upper.z);
70 return movemask((lower_x <= upper_x) & (lower_y <= upper_y) & (lower_z <= upper_z));
71 }
72
73 bool intersect_triangle_triangle (Scene* scene0, unsigned geomID0, unsigned primID0, Scene* scene1, unsigned geomID1, unsigned primID1)
74 {
75 CSTAT(bvh_collide_prim_intersections1++);
76 const TriangleMesh* mesh0 = scene0->get<TriangleMesh>(geomID0);
77 const TriangleMesh* mesh1 = scene1->get<TriangleMesh>(geomID1);
78 const TriangleMesh::Triangle& tri0 = mesh0->triangle(primID0);
79 const TriangleMesh::Triangle& tri1 = mesh1->triangle(primID1);
80
81 /* special culling for scene intersection with itself */
82 if (scene0 == scene1 && geomID0 == geomID1)
83 {
84 /* ignore self intersections */
85 if (primID0 == primID1)
86 return false;
87 }
88 CSTAT(bvh_collide_prim_intersections2++);
89
90 if (scene0 == scene1 && geomID0 == geomID1)
91 {
92 /* ignore intersection with topological neighbors */
93 const vint4 t0(tri0.v[0],tri0.v[1],tri0.v[2],tri0.v[2]);
94 if (any(vint4(tri1.v[0]) == t0)) return false;
95 if (any(vint4(tri1.v[1]) == t0)) return false;
96 if (any(vint4(tri1.v[2]) == t0)) return false;
97 }
98 CSTAT(bvh_collide_prim_intersections3++);
99
100 const Vec3fa a0 = mesh0->vertex(tri0.v[0]);
101 const Vec3fa a1 = mesh0->vertex(tri0.v[1]);
102 const Vec3fa a2 = mesh0->vertex(tri0.v[2]);
103 const Vec3fa b0 = mesh1->vertex(tri1.v[0]);
104 const Vec3fa b1 = mesh1->vertex(tri1.v[1]);
105 const Vec3fa b2 = mesh1->vertex(tri1.v[2]);
106
107 return TriangleTriangleIntersector::intersect_triangle_triangle(a0,a1,a2,b0,b1,b2);
108 }
109
110 template<int N>
111 __forceinline void BVHNColliderUserGeom<N>::processLeaf(NodeRef node0, NodeRef node1)
112 {
113 Collision collisions[16];
114 size_t num_collisions = 0;
115
116 size_t N0; Object* leaf0 = (Object*) node0.leaf(N0);
117 size_t N1; Object* leaf1 = (Object*) node1.leaf(N1);
118 for (size_t i=0; i<N0; i++) {
119 for (size_t j=0; j<N1; j++) {
120 const unsigned geomID0 = leaf0[i].geomID();
121 const unsigned primID0 = leaf0[i].primID();
122 const unsigned geomID1 = leaf1[j].geomID();
123 const unsigned primID1 = leaf1[j].primID();
124 if (this->scene0 == this->scene1 && geomID0 == geomID1 && primID0 == primID1) continue;
125 collisions[num_collisions++] = Collision(geomID0,primID0,geomID1,primID1);
126 if (num_collisions == 16) {
127 this->callback(this->userPtr,(RTCCollision*)&collisions,num_collisions);
128 num_collisions = 0;
129 }
130 }
131 }
132 if (num_collisions)
133 this->callback(this->userPtr,(RTCCollision*)&collisions,num_collisions);
134 }
135
136 template<int N>
137 void BVHNCollider<N>::collide_recurse(NodeRef ref0, const BBox3fa& bounds0, NodeRef ref1, const BBox3fa& bounds1, size_t depth0, size_t depth1)
138 {
139 CSTAT(bvh_collide_traversal_steps++);
140 if (unlikely(ref0.isLeaf())) {
141 if (unlikely(ref1.isLeaf())) {
142 CSTAT(bvh_collide_leaf_pairs++);
143 processLeaf(ref0,ref1);
144 return;
145 } else goto recurse_node1;
146
147 } else {
148 if (unlikely(ref1.isLeaf())) {
149 goto recurse_node0;
150 } else {
151 if (area(bounds0) > area(bounds1)) {
152 goto recurse_node0;
153 }
154 else {
155 goto recurse_node1;
156 }
157 }
158 }
159
160 {
161 recurse_node0:
162 AABBNode* node0 = ref0.getAABBNode();
163 size_t mask = overlap<N>(bounds1,*node0);
164 //for (size_t m=mask, i=bsf(m); m!=0; m=btc(m,i), i=bsf(m)) {
165 //for (size_t i=0; i<N; i++) {
166#if 0
167 if (depth0 < parallel_depth_threshold)
168 {
169 parallel_for(size_t(N), [&] ( size_t i ) {
170 if (mask & ( 1 << i)) {
171 BVHN<N>::prefetch(node0->child(i),BVH_FLAG_ALIGNED_NODE);
172 collide_recurse(node0->child(i),node0->bounds(i),ref1,bounds1,depth0+1,depth1);
173 }
174 });
175 }
176 else
177#endif
178 {
179 for (size_t m=mask, i=bsf(m); m!=0; m=btc(m,i), i=bsf(m)) {
180 BVHN<N>::prefetch(node0->child(i),BVH_FLAG_ALIGNED_NODE);
181 collide_recurse(node0->child(i),node0->bounds(i),ref1,bounds1,depth0+1,depth1);
182 }
183 }
184 return;
185 }
186
187 {
188 recurse_node1:
189 AABBNode* node1 = ref1.getAABBNode();
190 size_t mask = overlap<N>(bounds0,*node1);
191 //for (size_t m=mask, i=bsf(m); m!=0; m=btc(m,i), i=bsf(m)) {
192 //for (size_t i=0; i<N; i++) {
193#if 0
194 if (depth1 < parallel_depth_threshold)
195 {
196 parallel_for(size_t(N), [&] ( size_t i ) {
197 if (mask & ( 1 << i)) {
198 BVHN<N>::prefetch(node1->child(i),BVH_FLAG_ALIGNED_NODE);
199 collide_recurse(ref0,bounds0,node1->child(i),node1->bounds(i),depth0,depth1+1);
200 }
201 });
202 }
203 else
204#endif
205 {
206 for (size_t m=mask, i=bsf(m); m!=0; m=btc(m,i), i=bsf(m)) {
207 BVHN<N>::prefetch(node1->child(i),BVH_FLAG_ALIGNED_NODE);
208 collide_recurse(ref0,bounds0,node1->child(i),node1->bounds(i),depth0,depth1+1);
209 }
210 }
211 return;
212 }
213 }
214
215 template<int N>
216 void BVHNCollider<N>::split(const CollideJob& job, jobvector& jobs)
217 {
218 if (unlikely(job.ref0.isLeaf())) {
219 if (unlikely(job.ref1.isLeaf())) {
220 jobs.push_back(job);
221 return;
222 } else goto recurse_node1;
223 } else {
224 if (unlikely(job.ref1.isLeaf())) {
225 goto recurse_node0;
226 } else {
227 if (area(job.bounds0) > area(job.bounds1)) {
228 goto recurse_node0;
229 }
230 else {
231 goto recurse_node1;
232 }
233 }
234 }
235
236 {
237 recurse_node0:
238 const AABBNode* node0 = job.ref0.getAABBNode();
239 size_t mask = overlap<N>(job.bounds1,*node0);
240 for (size_t m=mask, i=bsf(m); m!=0; m=btc(m,i), i=bsf(m)) {
241 jobs.push_back(CollideJob(node0->child(i),node0->bounds(i),job.depth0+1,job.ref1,job.bounds1,job.depth1));
242 }
243 return;
244 }
245
246 {
247 recurse_node1:
248 const AABBNode* node1 = job.ref1.getAABBNode();
249 size_t mask = overlap<N>(job.bounds0,*node1);
250 for (size_t m=mask, i=bsf(m); m!=0; m=btc(m,i), i=bsf(m)) {
251 jobs.push_back(CollideJob(job.ref0,job.bounds0,job.depth0,node1->child(i),node1->bounds(i),job.depth1+1));
252 }
253 return;
254 }
255 }
256
257 template<int N>
258 void BVHNCollider<N>::collide_recurse_entry(NodeRef ref0, const BBox3fa& bounds0, NodeRef ref1, const BBox3fa& bounds1)
259 {
260 CSTAT(bvh_collide_traversal_steps = 0);
261 CSTAT(bvh_collide_leaf_pairs = 0);
262 CSTAT(bvh_collide_leaf_iterations = 0);
263 CSTAT(bvh_collide_prim_intersections1 = 0);
264 CSTAT(bvh_collide_prim_intersections2 = 0);
265 CSTAT(bvh_collide_prim_intersections3 = 0);
266 CSTAT(bvh_collide_prim_intersections4 = 0);
267 CSTAT(bvh_collide_prim_intersections5 = 0);
268 CSTAT(bvh_collide_prim_intersections = 0);
269#if 0
270 collide_recurse(ref0,bounds0,ref1,bounds1,0,0);
271#else
272 const int M = 2048;
273 jobvector jobs[2];
274 jobs[0].reserve(M);
275 jobs[1].reserve(M);
276 jobs[0].push_back(CollideJob(ref0,bounds0,0,ref1,bounds1,0));
277 int source = 0;
278 int target = 1;
279
280 /* try to split job until job list is full */
281 while (jobs[source].size()+8 <= M)
282 {
283 for (size_t i=0; i<jobs[source].size(); i++)
284 {
285 const CollideJob& job = jobs[source][i];
286 size_t remaining = jobs[source].size()-i;
287 if (jobs[target].size()+remaining+8 > M) {
288 jobs[target].push_back(job);
289 } else {
290 split(job,jobs[target]);
291 }
292 }
293
294 /* stop splitting jobs if we reached only leaves and cannot make progress anymore */
295 if (jobs[target].size() == jobs[source].size())
296 break;
297
298 jobs[source].resize(0);
299 std::swap(source,target);
300 }
301
302 /* parallel processing of all jobs */
303 parallel_for(size_t(jobs[source].size()), [&] ( size_t i ) {
304 CollideJob& j = jobs[source][i];
305 collide_recurse(j.ref0,j.bounds0,j.ref1,j.bounds1,j.depth0,j.depth1);
306 });
307
308
309#endif
310 CSTAT(PRINT(bvh_collide_traversal_steps));
311 CSTAT(PRINT(bvh_collide_leaf_pairs));
312 CSTAT(PRINT(bvh_collide_leaf_iterations));
313 CSTAT(PRINT(bvh_collide_prim_intersections1));
314 CSTAT(PRINT(bvh_collide_prim_intersections2));
315 CSTAT(PRINT(bvh_collide_prim_intersections3));
316 CSTAT(PRINT(bvh_collide_prim_intersections4));
317 CSTAT(PRINT(bvh_collide_prim_intersections5));
318 CSTAT(PRINT(bvh_collide_prim_intersections));
319 }
320
321 template<int N>
322 void BVHNColliderUserGeom<N>::collide(BVH* __restrict__ bvh0, BVH* __restrict__ bvh1, RTCCollideFunc callback, void* userPtr)
323 {
324 BVHNColliderUserGeom<N>(bvh0->scene,bvh1->scene,callback,userPtr).
325 collide_recurse_entry(bvh0->root,bvh0->bounds.bounds(),bvh1->root,bvh1->bounds.bounds());
326 }
327
328#if defined (EMBREE_LOWEST_ISA)
329 struct collision_regression_test : public RegressionTest
330 {
331 collision_regression_test(const char* name) : RegressionTest(name) {
332 registerRegressionTest(this);
333 }
334
335 bool run ()
336 {
337 bool passed = true;
338 passed &= TriangleTriangleIntersector::intersect_triangle_triangle (Vec3fa(-0.008815f, 0.041848f, -2.49875e-06f), Vec3fa(-0.008276f, 0.053318f, -2.49875e-06f), Vec3fa(0.003023f, 0.048969f, -2.49875e-06f),
339 Vec3fa(0.00245f, 0.037612f, -2.49875e-06f), Vec3fa(0.01434f, 0.042634f, -2.49875e-06f), Vec3fa(0.013499f, 0.031309f, -2.49875e-06f)) == false;
340 passed &= TriangleTriangleIntersector::intersect_triangle_triangle (Vec3fa(0,0,0),Vec3fa(1,0,0),Vec3fa(0,1,0), Vec3fa(0,0,0),Vec3fa(1,0,0),Vec3fa(0,1,0)) == true;
341 passed &= TriangleTriangleIntersector::intersect_triangle_triangle (Vec3fa(0,0,0),Vec3fa(1,0,0),Vec3fa(0,1,0), Vec3fa(0,0,1),Vec3fa(1,0,1),Vec3fa(0,1,1)) == false;
342 passed &= TriangleTriangleIntersector::intersect_triangle_triangle (Vec3fa(0,0,0),Vec3fa(1,0,0),Vec3fa(0,1,0), Vec3fa(0,0,1),Vec3fa(1,0,0),Vec3fa(0,1,0)) == true;
343 passed &= TriangleTriangleIntersector::intersect_triangle_triangle (Vec3fa(0,0,0),Vec3fa(1,0,0),Vec3fa(0,1,0), Vec3fa(0,0,0),Vec3fa(1,0,1),Vec3fa(0,1,1)) == true;
344 passed &= TriangleTriangleIntersector::intersect_triangle_triangle (Vec3fa(0,0,0),Vec3fa(1,0,0),Vec3fa(0,1,0), Vec3fa(0.1f,0.1f,0),Vec3fa(1,0,1),Vec3fa(0,1,1)) == true;
345 passed &= TriangleTriangleIntersector::intersect_triangle_triangle (Vec3fa(0,0,0),Vec3fa(1,0,0),Vec3fa(0,1,0), Vec3fa(0.1f,0.1f,-0.1f),Vec3fa(1,0,1),Vec3fa(0,1,1)) == true;
346 passed &= TriangleTriangleIntersector::intersect_triangle_triangle (Vec3fa(0,0,0),Vec3fa(1,0,0),Vec3fa(0,1,0), Vec3fa(0,0,0),Vec3fa(1,0,0),Vec3fa(0,1,0)) == true;
347 passed &= TriangleTriangleIntersector::intersect_triangle_triangle (Vec3fa(0,0,0),Vec3fa(1,0,0),Vec3fa(0,1,0), Vec3fa(0,0,0),Vec3fa(0.5f,0,0),Vec3fa(0,0.5f,0)) == true;
348 passed &= TriangleTriangleIntersector::intersect_triangle_triangle (Vec3fa(0,0,0),Vec3fa(1,0,0),Vec3fa(0,1,0), Vec3fa(0.1f,0.1f,0),Vec3fa(0.5f,0,0),Vec3fa(0,0.5f,0)) == true;
349 passed &= TriangleTriangleIntersector::intersect_triangle_triangle (Vec3fa(0,0,0),Vec3fa(1,0,0),Vec3fa(0,1,0), Vec3fa(0.1f,0.1f,0),Vec3fa(0.5f,0.1f,0),Vec3fa(0.1f,0.5f,0)) == true;
350 passed &= TriangleTriangleIntersector::intersect_triangle_triangle (Vec3fa(0,0,0),Vec3fa(1,0,0),Vec3fa(0,1,0), Vec3fa(0.1f,-0.1f,0),Vec3fa(0.5f,0.1f,0),Vec3fa(0.1f,0.5f,0)) == true;
351 passed &= TriangleTriangleIntersector::intersect_triangle_triangle (Vec3fa(0,0,0),Vec3fa(1,0,0),Vec3fa(0,1,0), Vec3fa(-0.1f,0.1f,0),Vec3fa(0.5f,0.1f,0),Vec3fa(0.1f,0.5f,0)) == true;
352 passed &= TriangleTriangleIntersector::intersect_triangle_triangle (Vec3fa(0,0,0),Vec3fa(1,0,0),Vec3fa(0,1,0),
353 Vec3fa(-1,1,0) + Vec3fa(0,0,0),Vec3fa(-1,1,0) + Vec3fa(0.1f,0,0),Vec3fa(-1,1,0) + Vec3fa(0,0.1f,0)) == false;
354 passed &= TriangleTriangleIntersector::intersect_triangle_triangle (Vec3fa(0,0,0),Vec3fa(1,0,0),Vec3fa(0,1,0),
355 Vec3fa( 2,0.5f,0) + Vec3fa(0,0,0),Vec3fa( 2,0.5f,0) + Vec3fa(0.1f,0,0),Vec3fa( 2,0.5f,0) + Vec3fa(0,0.1f,0)) == false;
356 passed &= TriangleTriangleIntersector::intersect_triangle_triangle (Vec3fa(0,0,0),Vec3fa(1,0,0),Vec3fa(0,1,0),
357 Vec3fa(0.5f,-2.0f,0) + Vec3fa(0,0,0),Vec3fa(0.5f,-2.0f,0) + Vec3fa(0.1f,0,0),Vec3fa(0.5f,-2.0f,0) + Vec3fa(0,0.1f,0)) == false;
358 return passed;
359 }
360 };
361
362 collision_regression_test collision_regression("collision_regression_test");
363#endif
364
365 ////////////////////////////////////////////////////////////////////////////////
366 /// Collider Definitions
367 ////////////////////////////////////////////////////////////////////////////////
368
369 DEFINE_COLLIDER(BVH4ColliderUserGeom,BVHNColliderUserGeom<4>);
370
371#if defined(__AVX__)
372 DEFINE_COLLIDER(BVH8ColliderUserGeom,BVHNColliderUserGeom<8>);
373#endif
374 }
375}
376