1 //============================================================================
2 // Copyright (c) Kitware, Inc.
3 // All rights reserved.
4 // See LICENSE.txt for details.
5 // This software is distributed WITHOUT ANY WARRANTY; without even
6 // the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
7 // PURPOSE. See the above copyright notice for more information.
8 //
9 // Copyright 2015 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
10 // Copyright 2015 UT-Battelle, LLC.
11 // Copyright 2015 Los Alamos National Security.
12 //
13 // Under the terms of Contract DE-NA0003525 with NTESS,
14 // the U.S. Government retains certain rights in this software.
15 //
16 // Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
17 // Laboratory (LANL), the U.S. Government retains certain rights in
18 // this software.
19 //============================================================================
20
21 #include <math.h>
22
23 #include <vtkm/Math.h>
24 #include <vtkm/VectorAnalysis.h>
25
26 #include <vtkm/cont/Algorithm.h>
27 #include <vtkm/cont/DeviceAdapter.h>
28 #include <vtkm/cont/DeviceAdapterAlgorithm.h>
29 #include <vtkm/cont/RuntimeDeviceTracker.h>
30 #include <vtkm/cont/Timer.h>
31 #include <vtkm/cont/TryExecute.h>
32
33 #include <vtkm/cont/AtomicArray.h>
34
35 #include <vtkm/rendering/raytracing/BoundingVolumeHierarchy.h>
36 #include <vtkm/rendering/raytracing/Logger.h>
37 #include <vtkm/rendering/raytracing/MortonCodes.h>
38 #include <vtkm/rendering/raytracing/RayTracingTypeDefs.h>
39 #include <vtkm/rendering/raytracing/Worklets.h>
40
41 #include <vtkm/worklet/DispatcherMapField.h>
42 #include <vtkm/worklet/WorkletMapField.h>
43
44 #define AABB_EPSILON 0.00001f
45 namespace vtkm
46 {
47 namespace rendering
48 {
49 namespace raytracing
50 {
51 namespace detail
52 {
53
54 class LinearBVHBuilder
55 {
56 public:
57 class CountingIterator;
58
59 template <typename Device>
60 class GatherFloat32;
61
62 template <typename Device>
63 class GatherVecCast;
64
65 class CreateLeafs;
66
67 class BVHData;
68
69 template <typename Device>
70 class PropagateAABBs;
71
72 template <typename Device>
73 class TreeBuilder;
74
75 VTKM_CONT
LinearBVHBuilder()76 LinearBVHBuilder() {}
77
78 template <typename Device>
79 VTKM_CONT void SortAABBS(BVHData& bvh, Device vtkmNotUsed(device), bool);
80
81 template <typename Device>
82 VTKM_CONT void BuildHierarchy(BVHData& bvh);
83
84 template <typename Device>
85 VTKM_CONT void RunOnDevice(LinearBVH& linearBVH, Device device);
86 }; // class LinearBVHBuilder
87
88 class LinearBVHBuilder::CountingIterator : public vtkm::worklet::WorkletMapField
89 {
90 public:
91 VTKM_CONT
CountingIterator()92 CountingIterator() {}
93 using ControlSignature = void(FieldOut<>);
94 using ExecutionSignature = void(WorkIndex, _1);
95 VTKM_EXEC
operator ()(const vtkm::Id & index,vtkm::Id & outId) const96 void operator()(const vtkm::Id& index, vtkm::Id& outId) const { outId = index; }
97 }; //class countingIterator
98
99 template <typename Device>
100 class LinearBVHBuilder::GatherFloat32 : public vtkm::worklet::WorkletMapField
101 {
102 private:
103 using FloatArrayHandle = typename vtkm::cont::ArrayHandle<vtkm::Float32>;
104 using PortalConst = typename FloatArrayHandle::ExecutionTypes<Device>::PortalConst;
105 using Portal = typename FloatArrayHandle::ExecutionTypes<Device>::Portal;
106 PortalConst InputPortal;
107 Portal OutputPortal;
108
109 public:
110 VTKM_CONT
GatherFloat32(const FloatArrayHandle & inputPortal,FloatArrayHandle & outputPortal,const vtkm::Id & size)111 GatherFloat32(const FloatArrayHandle& inputPortal,
112 FloatArrayHandle& outputPortal,
113 const vtkm::Id& size)
114 : InputPortal(inputPortal.PrepareForInput(Device()))
115 {
116 this->OutputPortal = outputPortal.PrepareForOutput(size, Device());
117 }
118 using ControlSignature = void(FieldIn<>);
119 using ExecutionSignature = void(WorkIndex, _1);
120 VTKM_EXEC
operator ()(const vtkm::Id & outIndex,const vtkm::Id & inIndex) const121 void operator()(const vtkm::Id& outIndex, const vtkm::Id& inIndex) const
122 {
123 OutputPortal.Set(outIndex, InputPortal.Get(inIndex));
124 }
125 }; //class GatherFloat
126
127 class LinearBVHBuilder::CreateLeafs : public vtkm::worklet::WorkletMapField
128 {
129
130 public:
131 VTKM_CONT
CreateLeafs()132 CreateLeafs() {}
133
134 typedef void ControlSignature(FieldIn<>, WholeArrayOut<>);
135 typedef void ExecutionSignature(_1, _2, WorkIndex);
136
137 template <typename LeafPortalType>
operator ()(const vtkm::Id & dataIndex,LeafPortalType & leafs,const vtkm::Id & index) const138 VTKM_EXEC void operator()(const vtkm::Id& dataIndex,
139 LeafPortalType& leafs,
140 const vtkm::Id& index) const
141 {
142 const vtkm::Id offset = index * 2;
143 leafs.Set(offset, 1); // number of primitives
144 leafs.Set(offset + 1, dataIndex); // number of primitives
145 }
146 }; //class createLeafs
147
148 template <typename Device>
149 class LinearBVHBuilder::GatherVecCast : public vtkm::worklet::WorkletMapField
150 {
151 private:
152 using Vec4IdArrayHandle = typename vtkm::cont::ArrayHandle<vtkm::Vec<vtkm::Id, 4>>;
153 using Vec4IntArrayHandle = typename vtkm::cont::ArrayHandle<vtkm::Vec<vtkm::Int32, 4>>;
154 using PortalConst = typename Vec4IdArrayHandle::ExecutionTypes<Device>::PortalConst;
155 using Portal = typename Vec4IntArrayHandle::ExecutionTypes<Device>::Portal;
156
157 private:
158 PortalConst InputPortal;
159 Portal OutputPortal;
160
161 public:
162 VTKM_CONT
GatherVecCast(const Vec4IdArrayHandle & inputPortal,Vec4IntArrayHandle & outputPortal,const vtkm::Id & size)163 GatherVecCast(const Vec4IdArrayHandle& inputPortal,
164 Vec4IntArrayHandle& outputPortal,
165 const vtkm::Id& size)
166 : InputPortal(inputPortal.PrepareForInput(Device()))
167 {
168 this->OutputPortal = outputPortal.PrepareForOutput(size, Device());
169 }
170 using ControlSignature = void(FieldIn<>);
171 using ExecutionSignature = void(WorkIndex, _1);
172 VTKM_EXEC
operator ()(const vtkm::Id & outIndex,const vtkm::Id & inIndex) const173 void operator()(const vtkm::Id& outIndex, const vtkm::Id& inIndex) const
174 {
175 OutputPortal.Set(outIndex, InputPortal.Get(inIndex));
176 }
177 }; //class GatherVec3Id
178
179 class LinearBVHBuilder::BVHData
180 {
181 public:
182 vtkm::cont::ArrayHandle<vtkm::UInt32> mortonCodes;
183 vtkm::cont::ArrayHandle<vtkm::Id> parent;
184 vtkm::cont::ArrayHandle<vtkm::Id> leftChild;
185 vtkm::cont::ArrayHandle<vtkm::Id> rightChild;
186 vtkm::cont::ArrayHandle<vtkm::Id> leafs;
187 vtkm::cont::ArrayHandle<vtkm::Bounds> innerBounds;
188 vtkm::cont::ArrayHandleCounting<vtkm::Id> leafOffsets;
189 AABBs& AABB;
190
191 template <typename Device>
BVHData(vtkm::Id numPrimitives,AABBs & aabbs,Device vtkmNotUsed (device))192 VTKM_CONT BVHData(vtkm::Id numPrimitives, AABBs& aabbs, Device vtkmNotUsed(device))
193 : leafOffsets(0, 2, numPrimitives)
194 , AABB(aabbs)
195 , NumPrimitives(numPrimitives)
196 {
197 InnerNodeCount = NumPrimitives - 1;
198 vtkm::Id size = NumPrimitives + InnerNodeCount;
199
200 parent.PrepareForOutput(size, Device());
201 leftChild.PrepareForOutput(InnerNodeCount, Device());
202 rightChild.PrepareForOutput(InnerNodeCount, Device());
203 innerBounds.PrepareForOutput(InnerNodeCount, Device());
204 mortonCodes.PrepareForOutput(NumPrimitives, Device());
205 }
206
207 VTKM_CONT
~BVHData()208 ~BVHData() {}
209
210 VTKM_CONT
GetNumberOfPrimitives() const211 vtkm::Id GetNumberOfPrimitives() const { return NumPrimitives; }
212 VTKM_CONT
GetNumberOfInnerNodes() const213 vtkm::Id GetNumberOfInnerNodes() const { return InnerNodeCount; }
214
215 private:
216 vtkm::Id NumPrimitives;
217 vtkm::Id InnerNodeCount;
218
219 }; // class BVH
220
221 template <typename Device>
222 class LinearBVHBuilder::PropagateAABBs : public vtkm::worklet::WorkletMapField
223 {
224 private:
225 using IdArrayHandle = typename vtkm::cont::ArrayHandle<vtkm::Id>;
226 using Int8Handle = typename vtkm::cont::ArrayHandle<vtkm::Int8>;
227 using Float2ArrayHandle = typename vtkm::cont::ArrayHandle<Vec<vtkm::Float32, 2>>;
228 using VecInt2Handle = typename vtkm::cont::ArrayHandle<vtkm::Vec<vtkm::Int32, 2>>;
229 using Float4ArrayHandle = typename vtkm::cont::ArrayHandle<Vec<vtkm::Float32, 4>>;
230
231 using IdConstPortal = typename IdArrayHandle::ExecutionTypes<Device>::PortalConst;
232 using Float2ArrayPortal = typename Float2ArrayHandle::ExecutionTypes<Device>::Portal;
233 using Int2ArrayPortal = typename VecInt2Handle::ExecutionTypes<Device>::Portal;
234 using Int8ArrayPortal = typename Int8Handle::ExecutionTypes<Device>::Portal;
235 using Float4ArrayPortal = typename Float4ArrayHandle::ExecutionTypes<Device>::Portal;
236
237 Float4ArrayPortal FlatBVH;
238 IdConstPortal Parents;
239 IdConstPortal LeftChildren;
240 IdConstPortal RightChildren;
241 vtkm::Int32 LeafCount;
242 vtkm::exec::AtomicArrayExecutionObject<vtkm::Int32, Device> Counters;
243
244 public:
245 VTKM_CONT
PropagateAABBs(IdArrayHandle & parents,IdArrayHandle & leftChildren,IdArrayHandle & rightChildren,vtkm::Int32 leafCount,Float4ArrayHandle flatBVH,const vtkm::cont::AtomicArray<vtkm::Int32> & counters)246 PropagateAABBs(IdArrayHandle& parents,
247 IdArrayHandle& leftChildren,
248 IdArrayHandle& rightChildren,
249 vtkm::Int32 leafCount,
250 Float4ArrayHandle flatBVH,
251 const vtkm::cont::AtomicArray<vtkm::Int32>& counters)
252 : Parents(parents.PrepareForInput(Device()))
253 , LeftChildren(leftChildren.PrepareForInput(Device()))
254 , RightChildren(rightChildren.PrepareForInput(Device()))
255 , LeafCount(leafCount)
256 , Counters(counters.PrepareForExecution(Device()))
257
258 {
259 this->FlatBVH = flatBVH.PrepareForOutput((LeafCount - 1) * 4, Device());
260 }
261 using ControlSignature = void(WholeArrayIn<Scalar>,
262 WholeArrayIn<Scalar>,
263 WholeArrayIn<Scalar>,
264 WholeArrayIn<Scalar>,
265 WholeArrayIn<Scalar>,
266 WholeArrayIn<Scalar>,
267 WholeArrayIn<>);
268 using ExecutionSignature = void(WorkIndex, _1, _2, _3, _4, _5, _6, _7);
269
270 template <typename InputPortalType, typename OffsetPortalType>
operator ()(const vtkm::Id workIndex,const InputPortalType & xmin,const InputPortalType & ymin,const InputPortalType & zmin,const InputPortalType & xmax,const InputPortalType & ymax,const InputPortalType & zmax,const OffsetPortalType & leafOffsets) const271 VTKM_EXEC_CONT void operator()(const vtkm::Id workIndex,
272 const InputPortalType& xmin,
273 const InputPortalType& ymin,
274 const InputPortalType& zmin,
275 const InputPortalType& xmax,
276 const InputPortalType& ymax,
277 const InputPortalType& zmax,
278 const OffsetPortalType& leafOffsets) const
279 {
280 //move up into the inner nodes
281 vtkm::Id currentNode = LeafCount - 1 + workIndex;
282 vtkm::Vec<vtkm::Id, 2> childVector;
283 while (currentNode != 0)
284 {
285 currentNode = Parents.Get(currentNode);
286
287 vtkm::Int32 oldCount = Counters.Add(currentNode, 1);
288 if (oldCount == 0)
289 return;
290 vtkm::Id currentNodeOffset = currentNode * 4;
291 childVector[0] = LeftChildren.Get(currentNode);
292 childVector[1] = RightChildren.Get(currentNode);
293 if (childVector[0] > (LeafCount - 2))
294 {
295 //our left child is a leaf, so just grab the AABB
296 //and set it in the current node
297 childVector[0] = childVector[0] - LeafCount + 1;
298
299 vtkm::Vec<vtkm::Float32, 4>
300 first4Vec; // = FlatBVH.Get(currentNode); only this one needs effects this
301
302 first4Vec[0] = xmin.Get(childVector[0]);
303 first4Vec[1] = ymin.Get(childVector[0]);
304 first4Vec[2] = zmin.Get(childVector[0]);
305 first4Vec[3] = xmax.Get(childVector[0]);
306 FlatBVH.Set(currentNodeOffset, first4Vec);
307
308 vtkm::Vec<vtkm::Float32, 4> second4Vec = FlatBVH.Get(currentNodeOffset + 1);
309 second4Vec[0] = ymax.Get(childVector[0]);
310 second4Vec[1] = zmax.Get(childVector[0]);
311 FlatBVH.Set(currentNodeOffset + 1, second4Vec);
312 // set index to leaf
313 vtkm::Id leafIndex = leafOffsets.Get(childVector[0]);
314 childVector[0] = -(leafIndex + 1);
315 }
316 else
317 {
318 //our left child is an inner node, so gather
319 //both AABBs in the child and join them for
320 //the current node left AABB.
321 vtkm::Id child = childVector[0] * 4;
322
323 vtkm::Vec<vtkm::Float32, 4> cFirst4Vec = FlatBVH.Get(child);
324 vtkm::Vec<vtkm::Float32, 4> cSecond4Vec = FlatBVH.Get(child + 1);
325 vtkm::Vec<vtkm::Float32, 4> cThird4Vec = FlatBVH.Get(child + 2);
326
327 cFirst4Vec[0] = vtkm::Min(cFirst4Vec[0], cSecond4Vec[2]);
328 cFirst4Vec[1] = vtkm::Min(cFirst4Vec[1], cSecond4Vec[3]);
329 cFirst4Vec[2] = vtkm::Min(cFirst4Vec[2], cThird4Vec[0]);
330 cFirst4Vec[3] = vtkm::Max(cFirst4Vec[3], cThird4Vec[1]);
331 FlatBVH.Set(currentNodeOffset, cFirst4Vec);
332
333 vtkm::Vec<vtkm::Float32, 4> second4Vec = FlatBVH.Get(currentNodeOffset + 1);
334 second4Vec[0] = vtkm::Max(cSecond4Vec[0], cThird4Vec[2]);
335 second4Vec[1] = vtkm::Max(cSecond4Vec[1], cThird4Vec[3]);
336
337 FlatBVH.Set(currentNodeOffset + 1, second4Vec);
338 }
339
340 if (childVector[1] > (LeafCount - 2))
341 {
342 //our right child is a leaf, so just grab the AABB
343 //and set it in the current node
344 childVector[1] = childVector[1] - LeafCount + 1;
345
346
347 vtkm::Vec<vtkm::Float32, 4> second4Vec = FlatBVH.Get(currentNodeOffset + 1);
348
349 second4Vec[2] = xmin.Get(childVector[1]);
350 second4Vec[3] = ymin.Get(childVector[1]);
351 FlatBVH.Set(currentNodeOffset + 1, second4Vec);
352
353 vtkm::Vec<vtkm::Float32, 4> third4Vec;
354 third4Vec[0] = zmin.Get(childVector[1]);
355 third4Vec[1] = xmax.Get(childVector[1]);
356 third4Vec[2] = ymax.Get(childVector[1]);
357 third4Vec[3] = zmax.Get(childVector[1]);
358 FlatBVH.Set(currentNodeOffset + 2, third4Vec);
359
360 // set index to leaf
361 vtkm::Id leafIndex = leafOffsets.Get(childVector[1]);
362 childVector[1] = -(leafIndex + 1);
363 }
364 else
365 {
366 //our left child is an inner node, so gather
367 //both AABBs in the child and join them for
368 //the current node left AABB.
369 vtkm::Id child = childVector[1] * 4;
370
371 vtkm::Vec<vtkm::Float32, 4> cFirst4Vec = FlatBVH.Get(child);
372 vtkm::Vec<vtkm::Float32, 4> cSecond4Vec = FlatBVH.Get(child + 1);
373 vtkm::Vec<vtkm::Float32, 4> cThird4Vec = FlatBVH.Get(child + 2);
374
375 vtkm::Vec<vtkm::Float32, 4> second4Vec = FlatBVH.Get(currentNodeOffset + 1);
376 second4Vec[2] = vtkm::Min(cFirst4Vec[0], cSecond4Vec[2]);
377 second4Vec[3] = vtkm::Min(cFirst4Vec[1], cSecond4Vec[3]);
378 FlatBVH.Set(currentNodeOffset + 1, second4Vec);
379
380 cThird4Vec[0] = vtkm::Min(cFirst4Vec[2], cThird4Vec[0]);
381 cThird4Vec[1] = vtkm::Max(cFirst4Vec[3], cThird4Vec[1]);
382 cThird4Vec[2] = vtkm::Max(cSecond4Vec[0], cThird4Vec[2]);
383 cThird4Vec[3] = vtkm::Max(cSecond4Vec[1], cThird4Vec[3]);
384 FlatBVH.Set(currentNodeOffset + 2, cThird4Vec);
385 }
386 vtkm::Vec<vtkm::Float32, 4> fourth4Vec;
387 vtkm::Int32 leftChild =
388 static_cast<vtkm::Int32>((childVector[0] >= 0) ? childVector[0] * 4 : childVector[0]);
389 memcpy(&fourth4Vec[0], &leftChild, 4);
390 vtkm::Int32 rightChild =
391 static_cast<vtkm::Int32>((childVector[1] >= 0) ? childVector[1] * 4 : childVector[1]);
392 memcpy(&fourth4Vec[1], &rightChild, 4);
393 FlatBVH.Set(currentNodeOffset + 3, fourth4Vec);
394 }
395 }
396 }; //class PropagateAABBs
397
398 template <typename Device>
399 class LinearBVHBuilder::TreeBuilder : public vtkm::worklet::WorkletMapField
400 {
401 public:
402 using UIntArrayHandle = typename vtkm::cont::ArrayHandle<vtkm::UInt32>;
403 using IdArrayHandle = typename vtkm::cont::ArrayHandle<vtkm::Id>;
404 using UIntPortalType = typename UIntArrayHandle::ExecutionTypes<Device>::PortalConst;
405 using IdPortalType = typename IdArrayHandle::ExecutionTypes<Device>::Portal;
406
407 private:
408 UIntPortalType MortonCodePortal;
409 IdPortalType ParentPortal;
410 vtkm::Id LeafCount;
411 vtkm::Id InnerCount;
412 //TODO: get intrinsic support
413 VTKM_EXEC
CountLeadingZeros(vtkm::UInt32 & x) const414 inline vtkm::Int32 CountLeadingZeros(vtkm::UInt32& x) const
415 {
416 vtkm::UInt32 y;
417 vtkm::UInt32 n = 32;
418 y = x >> 16;
419 if (y != 0)
420 {
421 n = n - 16;
422 x = y;
423 }
424 y = x >> 8;
425 if (y != 0)
426 {
427 n = n - 8;
428 x = y;
429 }
430 y = x >> 4;
431 if (y != 0)
432 {
433 n = n - 4;
434 x = y;
435 }
436 y = x >> 2;
437 if (y != 0)
438 {
439 n = n - 2;
440 x = y;
441 }
442 y = x >> 1;
443 if (y != 0)
444 return vtkm::Int32(n - 2);
445 return vtkm::Int32(n - x);
446 }
447
448 // returns the count of largest shared prefix between
449 // two morton codes. Ties are broken by the indexes
450 // a and b.
451 //
452 // returns count of the largest binary prefix
453
454 VTKM_EXEC
delta(const vtkm::Int32 & a,const vtkm::Int32 & b) const455 inline vtkm::Int32 delta(const vtkm::Int32& a, const vtkm::Int32& b) const
456 {
457 bool tie = false;
458 bool outOfRange = (b < 0 || b > LeafCount - 1);
459 //still make the call but with a valid adderss
460 vtkm::Int32 bb = (outOfRange) ? 0 : b;
461 vtkm::UInt32 aCode = MortonCodePortal.Get(a);
462 vtkm::UInt32 bCode = MortonCodePortal.Get(bb);
463 //use xor to find where they differ
464 vtkm::UInt32 exOr = aCode ^ bCode;
465 tie = (exOr == 0);
466 //break the tie, a and b must always differ
467 exOr = tie ? vtkm::UInt32(a) ^ vtkm::UInt32(bb) : exOr;
468 vtkm::Int32 count = CountLeadingZeros(exOr);
469 if (tie)
470 count += 32;
471 count = (outOfRange) ? -1 : count;
472 return count;
473 }
474
475 public:
476 VTKM_CONT
TreeBuilder(const UIntArrayHandle & mortonCodesHandle,IdArrayHandle & parentHandle,const vtkm::Id & leafCount)477 TreeBuilder(const UIntArrayHandle& mortonCodesHandle,
478 IdArrayHandle& parentHandle,
479 const vtkm::Id& leafCount)
480 : MortonCodePortal(mortonCodesHandle.PrepareForInput(Device()))
481 , LeafCount(leafCount)
482 {
483 InnerCount = LeafCount - 1;
484 this->ParentPortal = parentHandle.PrepareForOutput(InnerCount + LeafCount, Device());
485 }
486 using ControlSignature = void(FieldOut<>, FieldOut<>);
487 using ExecutionSignature = void(WorkIndex, _1, _2);
488 VTKM_EXEC
operator ()(const vtkm::Id & index,vtkm::Id & leftChild,vtkm::Id & rightChild) const489 void operator()(const vtkm::Id& index, vtkm::Id& leftChild, vtkm::Id& rightChild) const
490 {
491 vtkm::Int32 idx = vtkm::Int32(index);
492 //something = MortonCodePortal.Get(index) + 1;
493 //determine range direction
494 vtkm::Int32 d = 0 > (delta(idx, idx + 1) - delta(idx, idx - 1)) ? -1 : 1;
495
496 //find upper bound for the length of the range
497 vtkm::Int32 minDelta = delta(idx, idx - d);
498 vtkm::Int32 lMax = 2;
499 while (delta(idx, idx + lMax * d) > minDelta)
500 lMax *= 2;
501
502 //binary search to find the lower bound
503 vtkm::Int32 l = 0;
504 for (int t = lMax / 2; t >= 1; t /= 2)
505 {
506 if (delta(idx, idx + (l + t) * d) > minDelta)
507 l += t;
508 }
509
510 vtkm::Int32 j = idx + l * d;
511 vtkm::Int32 deltaNode = delta(idx, j);
512 vtkm::Int32 s = 0;
513 vtkm::Float32 divFactor = 2.f;
514 //find the split position using a binary search
515 for (vtkm::Int32 t = (vtkm::Int32)ceil(vtkm::Float32(l) / divFactor);;
516 divFactor *= 2, t = (vtkm::Int32)ceil(vtkm::Float32(l) / divFactor))
517 {
518 if (delta(idx, idx + (s + t) * d) > deltaNode)
519 {
520 s += t;
521 }
522
523 if (t == 1)
524 break;
525 }
526
527 vtkm::Int32 split = idx + s * d + vtkm::Min(d, 0);
528 //assign parent/child pointers
529 if (vtkm::Min(idx, j) == split)
530 {
531 //leaf
532 ParentPortal.Set(split + InnerCount, idx);
533 leftChild = split + InnerCount;
534 }
535 else
536 {
537 //inner node
538 ParentPortal.Set(split, idx);
539 leftChild = split;
540 }
541
542
543 if (vtkm::Max(idx, j) == split + 1)
544 {
545 //leaf
546 ParentPortal.Set(split + InnerCount + 1, idx);
547 rightChild = split + InnerCount + 1;
548 }
549 else
550 {
551 ParentPortal.Set(split + 1, idx);
552 rightChild = split + 1;
553 }
554 }
555 }; // class TreeBuilder
556
557 template <typename Device>
SortAABBS(BVHData & bvh,Device vtkmNotUsed (device),bool singleAABB)558 VTKM_CONT void LinearBVHBuilder::SortAABBS(BVHData& bvh,
559 Device vtkmNotUsed(device),
560 bool singleAABB)
561 {
562
563 //create array of indexes to be sorted with morton codes
564 vtkm::cont::ArrayHandle<vtkm::Id> iterator;
565 iterator.PrepareForOutput(bvh.GetNumberOfPrimitives(), Device());
566 vtkm::worklet::DispatcherMapField<CountingIterator> iteratorDispatcher;
567 iteratorDispatcher.SetDevice(Device());
568 iteratorDispatcher.Invoke(iterator);
569
570 //sort the morton codes
571
572 vtkm::cont::DeviceAdapterAlgorithm<Device>::SortByKey(bvh.mortonCodes, iterator);
573
574 vtkm::Id arraySize = bvh.GetNumberOfPrimitives();
575 vtkm::cont::ArrayHandle<vtkm::Float32> temp1;
576 vtkm::cont::ArrayHandle<vtkm::Float32> temp2;
577
578
579 //tempStorage = new vtkm::cont::ArrayHandle<vtkm::Float32>();
580 //xmins
581 {
582 vtkm::worklet::DispatcherMapField<GatherFloat32<Device>> dispatcher(
583 GatherFloat32<Device>(bvh.AABB.xmins, temp1, arraySize));
584 dispatcher.SetDevice(Device());
585 dispatcher.Invoke(iterator);
586 }
587 temp2 = bvh.AABB.xmins;
588 bvh.AABB.xmins = temp1;
589 temp1 = temp2;
590
591 {
592 vtkm::worklet::DispatcherMapField<GatherFloat32<Device>> dispatcher(
593 GatherFloat32<Device>(bvh.AABB.ymins, temp1, arraySize));
594 dispatcher.SetDevice(Device());
595 dispatcher.Invoke(iterator);
596 }
597
598 temp2 = bvh.AABB.ymins;
599 bvh.AABB.ymins = temp1;
600 temp1 = temp2;
601 //zmins
602 {
603 vtkm::worklet::DispatcherMapField<GatherFloat32<Device>> dispatcher(
604 GatherFloat32<Device>(bvh.AABB.zmins, temp1, arraySize));
605 dispatcher.SetDevice(Device());
606 dispatcher.Invoke(iterator);
607 }
608
609 temp2 = bvh.AABB.zmins;
610 bvh.AABB.zmins = temp1;
611 temp1 = temp2;
612 //xmaxs
613 {
614 vtkm::worklet::DispatcherMapField<GatherFloat32<Device>> dispatcher(
615 GatherFloat32<Device>(bvh.AABB.xmaxs, temp1, arraySize));
616 dispatcher.SetDevice(Device());
617 dispatcher.Invoke(iterator);
618 }
619
620 temp2 = bvh.AABB.xmaxs;
621 bvh.AABB.xmaxs = temp1;
622 temp1 = temp2;
623 //ymaxs
624 {
625 vtkm::worklet::DispatcherMapField<GatherFloat32<Device>> dispatcher(
626 GatherFloat32<Device>(bvh.AABB.ymaxs, temp1, arraySize));
627 dispatcher.SetDevice(Device());
628 dispatcher.Invoke(iterator);
629 }
630
631 temp2 = bvh.AABB.ymaxs;
632 bvh.AABB.ymaxs = temp1;
633 temp1 = temp2;
634 //zmaxs
635 {
636 vtkm::worklet::DispatcherMapField<GatherFloat32<Device>> dispatcher(
637 GatherFloat32<Device>(bvh.AABB.zmaxs, temp1, arraySize));
638 dispatcher.SetDevice(Device());
639 dispatcher.Invoke(iterator);
640 }
641
642 temp2 = bvh.AABB.zmaxs;
643 bvh.AABB.zmaxs = temp1;
644 temp1 = temp2;
645
646 // Create the leaf references
647 bvh.leafs.PrepareForOutput(arraySize * 2, Device());
648 // we only actually have a single primitive, but the algorithm
649 // requires 2. Make sure they both point to the original
650 // primitive
651 if (singleAABB)
652 {
653 auto iterPortal = iterator.GetPortalControl();
654 for (int i = 0; i < 2; ++i)
655 {
656 iterPortal.Set(i, 0);
657 }
658 }
659 vtkm::worklet::DispatcherMapField<CreateLeafs> createDis;
660 createDis.SetDevice(Device());
661 createDis.Invoke(iterator, bvh.leafs);
662
663 } // method SortAABB
664
665 // Adding this as a template parameter to allow restricted types and
666 // storage for dynamic coordinate system to limit crazy code bloat and
667 // compile times.
668 //
669 template <typename Device>
RunOnDevice(LinearBVH & linearBVH,Device device)670 VTKM_CONT void LinearBVHBuilder::RunOnDevice(LinearBVH& linearBVH, Device device)
671 {
672 Logger* logger = Logger::GetInstance();
673 logger->OpenLogEntry("bvh_constuct");
674
675 vtkm::cont::Timer<Device> constructTimer;
676 //
677 //
678 // This algorithm needs at least 2 AABBs
679 //
680 bool singleAABB = false;
681 vtkm::Id numberOfAABBs = linearBVH.GetNumberOfAABBs();
682 if (numberOfAABBs == 1)
683 {
684 numberOfAABBs = 2;
685 singleAABB = true;
686 vtkm::Float32 xmin = linearBVH.AABB.xmins.GetPortalControl().Get(0);
687 vtkm::Float32 ymin = linearBVH.AABB.ymins.GetPortalControl().Get(0);
688 vtkm::Float32 zmin = linearBVH.AABB.zmins.GetPortalControl().Get(0);
689 vtkm::Float32 xmax = linearBVH.AABB.xmaxs.GetPortalControl().Get(0);
690 vtkm::Float32 ymax = linearBVH.AABB.ymaxs.GetPortalControl().Get(0);
691 vtkm::Float32 zmax = linearBVH.AABB.zmaxs.GetPortalControl().Get(0);
692
693 linearBVH.AABB.xmins.Allocate(2);
694 linearBVH.AABB.ymins.Allocate(2);
695 linearBVH.AABB.zmins.Allocate(2);
696 linearBVH.AABB.xmaxs.Allocate(2);
697 linearBVH.AABB.ymaxs.Allocate(2);
698 linearBVH.AABB.zmaxs.Allocate(2);
699 for (int i = 0; i < 2; ++i)
700 {
701 linearBVH.AABB.xmins.GetPortalControl().Set(i, xmin);
702 linearBVH.AABB.ymins.GetPortalControl().Set(i, ymin);
703 linearBVH.AABB.zmins.GetPortalControl().Set(i, zmin);
704 linearBVH.AABB.xmaxs.GetPortalControl().Set(i, xmax);
705 linearBVH.AABB.ymaxs.GetPortalControl().Set(i, ymax);
706 linearBVH.AABB.zmaxs.GetPortalControl().Set(i, zmax);
707 }
708 }
709
710
711 logger->AddLogData("bvh_num_aabbs", numberOfAABBs);
712
713 const vtkm::Id numBBoxes = numberOfAABBs;
714 BVHData bvh(numBBoxes, linearBVH.GetAABBs(), device);
715
716
717 vtkm::cont::Timer<Device> timer;
718 // Find the extent of all bounding boxes to generate normalization for morton codes
719 vtkm::Vec<vtkm::Float32, 3> minExtent(vtkm::Infinity32(), vtkm::Infinity32(), vtkm::Infinity32());
720 vtkm::Vec<vtkm::Float32, 3> maxExtent(
721 vtkm::NegativeInfinity32(), vtkm::NegativeInfinity32(), vtkm::NegativeInfinity32());
722 maxExtent[0] =
723 vtkm::cont::DeviceAdapterAlgorithm<Device>::Reduce(bvh.AABB.xmaxs, maxExtent[0], MaxValue());
724 maxExtent[1] =
725 vtkm::cont::DeviceAdapterAlgorithm<Device>::Reduce(bvh.AABB.ymaxs, maxExtent[1], MaxValue());
726 maxExtent[2] =
727 vtkm::cont::DeviceAdapterAlgorithm<Device>::Reduce(bvh.AABB.zmaxs, maxExtent[2], MaxValue());
728 minExtent[0] =
729 vtkm::cont::DeviceAdapterAlgorithm<Device>::Reduce(bvh.AABB.xmins, minExtent[0], MinValue());
730 minExtent[1] =
731 vtkm::cont::DeviceAdapterAlgorithm<Device>::Reduce(bvh.AABB.ymins, minExtent[1], MinValue());
732 minExtent[2] =
733 vtkm::cont::DeviceAdapterAlgorithm<Device>::Reduce(bvh.AABB.zmins, minExtent[2], MinValue());
734
735 linearBVH.TotalBounds.X.Min = minExtent[0];
736 linearBVH.TotalBounds.X.Max = maxExtent[0];
737 linearBVH.TotalBounds.Y.Min = minExtent[1];
738 linearBVH.TotalBounds.Y.Max = maxExtent[1];
739 linearBVH.TotalBounds.Z.Min = minExtent[2];
740 linearBVH.TotalBounds.Z.Max = maxExtent[2];
741
742 vtkm::Float64 time = timer.GetElapsedTime();
743 logger->AddLogData("calc_extents", time);
744 timer.Reset();
745
746 vtkm::Vec<vtkm::Float32, 3> deltaExtent = maxExtent - minExtent;
747 vtkm::Vec<vtkm::Float32, 3> inverseExtent;
748 for (int i = 0; i < 3; ++i)
749 {
750 inverseExtent[i] = (deltaExtent[i] == 0.f) ? 0 : 1.f / deltaExtent[i];
751 }
752
753 //Generate the morton codes
754 vtkm::worklet::DispatcherMapField<MortonCodeAABB> mortonDis(
755 MortonCodeAABB(inverseExtent, minExtent));
756 mortonDis.SetDevice(Device());
757 mortonDis.Invoke(bvh.AABB.xmins,
758 bvh.AABB.ymins,
759 bvh.AABB.zmins,
760 bvh.AABB.xmaxs,
761 bvh.AABB.ymaxs,
762 bvh.AABB.zmaxs,
763 bvh.mortonCodes);
764
765 time = timer.GetElapsedTime();
766 logger->AddLogData("morton_codes", time);
767 timer.Reset();
768
769 linearBVH.Allocate(bvh.GetNumberOfPrimitives(), Device());
770
771 SortAABBS(bvh, Device(), singleAABB);
772
773 time = timer.GetElapsedTime();
774 logger->AddLogData("sort_aabbs", time);
775 timer.Reset();
776
777 vtkm::worklet::DispatcherMapField<TreeBuilder<Device>> treeDis(
778 TreeBuilder<Device>(bvh.mortonCodes, bvh.parent, bvh.GetNumberOfPrimitives()));
779 treeDis.SetDevice(Device());
780 treeDis.Invoke(bvh.leftChild, bvh.rightChild);
781
782 time = timer.GetElapsedTime();
783 logger->AddLogData("build_tree", time);
784 timer.Reset();
785
786 const vtkm::Int32 primitiveCount = vtkm::Int32(bvh.GetNumberOfPrimitives());
787
788 vtkm::cont::ArrayHandle<vtkm::Int32> counters;
789 counters.PrepareForOutput(bvh.GetNumberOfPrimitives() - 1, Device());
790
791 vtkm::cont::ArrayHandleConstant<vtkm::Int32> zero(0, bvh.GetNumberOfPrimitives() - 1);
792 vtkm::cont::Algorithm::Copy(Device(), zero, counters);
793
794 vtkm::cont::AtomicArray<vtkm::Int32> atomicCounters(counters);
795
796
797 vtkm::worklet::DispatcherMapField<PropagateAABBs<Device>> propDis(PropagateAABBs<Device>(
798 bvh.parent, bvh.leftChild, bvh.rightChild, primitiveCount, linearBVH.FlatBVH, atomicCounters));
799 propDis.SetDevice(Device());
800 propDis.Invoke(bvh.AABB.xmins,
801 bvh.AABB.ymins,
802 bvh.AABB.zmins,
803 bvh.AABB.xmaxs,
804 bvh.AABB.ymaxs,
805 bvh.AABB.zmaxs,
806 bvh.leafOffsets);
807
808 linearBVH.Leafs = bvh.leafs;
809 time = timer.GetElapsedTime();
810 logger->AddLogData("propagate_aabbs", time);
811
812 time = constructTimer.GetElapsedTime();
813 logger->CloseLogEntry(time);
814 }
815 } //namespace detail
816
817 struct LinearBVH::ConstructFunctor
818 {
819 LinearBVH* Self;
820 VTKM_CONT
ConstructFunctorvtkm::rendering::raytracing::LinearBVH::ConstructFunctor821 ConstructFunctor(LinearBVH* self)
822 : Self(self)
823 {
824 }
825 template <typename Device>
operator ()vtkm::rendering::raytracing::LinearBVH::ConstructFunctor826 bool operator()(Device)
827 {
828 Self->ConstructOnDevice(Device());
829 return true;
830 }
831 };
832
LinearBVH()833 LinearBVH::LinearBVH()
834 : IsConstructed(false)
835 , CanConstruct(false){};
836
837 VTKM_CONT
LinearBVH(AABBs & aabbs)838 LinearBVH::LinearBVH(AABBs& aabbs)
839 : AABB(aabbs)
840 , IsConstructed(false)
841 , CanConstruct(true)
842 {
843 }
844
845 VTKM_CONT
LinearBVH(const LinearBVH & other)846 LinearBVH::LinearBVH(const LinearBVH& other)
847 : AABB(other.AABB)
848 , FlatBVH(other.FlatBVH)
849 , Leafs(other.Leafs)
850 , LeafCount(other.LeafCount)
851 , IsConstructed(other.IsConstructed)
852 , CanConstruct(other.CanConstruct)
853 {
854 }
855 template <typename Device>
Allocate(const vtkm::Id & leafCount,Device deviceAdapter)856 VTKM_CONT void LinearBVH::Allocate(const vtkm::Id& leafCount, Device deviceAdapter)
857 {
858 LeafCount = leafCount;
859 FlatBVH.PrepareForOutput((leafCount - 1) * 4, deviceAdapter);
860 }
861
Construct()862 void LinearBVH::Construct()
863 {
864 if (IsConstructed)
865 return;
866 if (!CanConstruct)
867 throw vtkm::cont::ErrorBadValue(
868 "Linear BVH: coordinates and triangles must be set before calling construct!");
869
870 ConstructFunctor functor(this);
871 vtkm::cont::TryExecute(functor);
872 IsConstructed = true;
873 }
874
875 VTKM_CONT
SetData(AABBs & aabbs)876 void LinearBVH::SetData(AABBs& aabbs)
877 {
878 AABB = aabbs;
879 IsConstructed = false;
880 CanConstruct = true;
881 }
882
883 template <typename Device>
ConstructOnDevice(Device device)884 void LinearBVH::ConstructOnDevice(Device device)
885 {
886 Logger* logger = Logger::GetInstance();
887 vtkm::cont::Timer<Device> timer;
888 logger->OpenLogEntry("bvh");
889 if (!CanConstruct)
890 throw vtkm::cont::ErrorBadValue(
891 "Linear BVH: coordinates and triangles must be set before calling construct!");
892 if (!IsConstructed)
893 {
894 detail::LinearBVHBuilder builder;
895 builder.RunOnDevice(*this, device);
896 IsConstructed = true;
897 }
898
899 vtkm::Float64 time = timer.GetElapsedTime();
900 logger->CloseLogEntry(time);
901 }
902
903 // explicitly export
904 template VTKM_RENDERING_EXPORT void LinearBVH::ConstructOnDevice<
905 vtkm::cont::DeviceAdapterTagSerial>(vtkm::cont::DeviceAdapterTagSerial);
906 #ifdef VTKM_ENABLE_TBB
907 template VTKM_RENDERING_EXPORT void LinearBVH::ConstructOnDevice<vtkm::cont::DeviceAdapterTagTBB>(
908 vtkm::cont::DeviceAdapterTagTBB);
909 #endif
910 #ifdef VTKM_ENABLE_OPENMP
911 template VTKM_CONT_EXPORT void LinearBVH::ConstructOnDevice<vtkm::cont::DeviceAdapterTagOpenMP>(
912 vtkm::cont::DeviceAdapterTagOpenMP);
913 #endif
914 #ifdef VTKM_ENABLE_CUDA
915 template VTKM_RENDERING_EXPORT void LinearBVH::ConstructOnDevice<vtkm::cont::DeviceAdapterTagCuda>(
916 vtkm::cont::DeviceAdapterTagCuda);
917 #endif
918
919 VTKM_CONT
GetIsConstructed() const920 bool LinearBVH::GetIsConstructed() const
921 {
922 return IsConstructed;
923 }
924
GetNumberOfAABBs() const925 vtkm::Id LinearBVH::GetNumberOfAABBs() const
926 {
927 return AABB.xmins.GetPortalConstControl().GetNumberOfValues();
928 }
929
GetAABBs()930 AABBs& LinearBVH::GetAABBs()
931 {
932 return AABB;
933 }
934 }
935 }
936 } // namespace vtkm::rendering::raytracing
937