OpenVDB 11.0.0
Loading...
Searching...
No Matches
GridStats.h
Go to the documentation of this file.
1// Copyright Contributors to the OpenVDB Project
2// SPDX-License-Identifier: MPL-2.0
3
4/*!
5 \file GridStats.h
6
7 \author Ken Museth
8
9 \date August 29, 2020
10
11 \brief Re-computes min/max/avg/var/bbox information for each node in a
12 pre-existing NanoVDB grid.
13*/
14
15#ifndef NANOVDB_GRIDSTATS_H_HAS_BEEN_INCLUDED
16#define NANOVDB_GRIDSTATS_H_HAS_BEEN_INCLUDED
17
18#include <nanovdb/NanoVDB.h>
19#include "Range.h"
20#include "ForEach.h"
21
22#ifdef NANOVDB_USE_TBB
23#include <tbb/parallel_reduce.h>
24#endif
25
26#if defined(__CUDACC__)
27#include <cuda/std/limits>// for cuda::std::numeric_limits
28#else
29#include <limits.h>// for std::numeric_limits
30#endif
31
32#include <atomic>
33#include <iostream>
34
35namespace nanovdb {
36
37/// @brief Grid flags which indicate what extra information is present in the grid buffer
38enum class StatsMode : uint32_t {
39 Disable = 0,// disable the computation of any type of statistics (obviously the FASTEST!)
40 BBox = 1,// only compute the bbox of active values per node and total activeVoxelCount
41 MinMax = 2,// additionally compute extrema values
42 All = 3,// compute all of the statics, i.e. bbox, min/max, average and standard deviation
43 Default = 3,// default computational mode for statistics
44 End = 4,
45};
46
47/// @brief Re-computes the min/max, stats and bbox information for an existing NanoVDB Grid
48///
49/// @param grid Grid whose stats to update
50/// @param mode Mode of computation for the statistics.
51template<typename BuildT>
53
54//================================================================================================
55
56template<typename ValueT, int Rank = TensorTraits<ValueT>::Rank>
57class Extrema;
58
59/// @brief Template specialization of Extrema on scalar value types, i.e. rank = 0
60template<typename ValueT>
61class Extrema<ValueT, 0>
62{
63protected:
64 ValueT mMin, mMax;
65
66public:
67 using ValueType = ValueT;
69#if defined(__CUDACC__)
70 : mMin(cuda::std::numeric_limits<ValueT>::max())
71 , mMax(cuda::std::numeric_limits<ValueT>::lowest())
72#else
73 : mMin(std::numeric_limits<ValueT>::max())
74 , mMax(std::numeric_limits<ValueT>::lowest())
75#endif
76 {
77 }
78 __hostdev__ Extrema(const ValueT& v)
79 : mMin(v)
80 , mMax(v)
81 {
82 }
83 __hostdev__ Extrema(const ValueT& a, const ValueT& b)
84 : mMin(a)
85 , mMax(b)
86 {
87 }
88 __hostdev__ Extrema& min(const ValueT& v)
89 {
90 if (v < mMin) mMin = v;
91 return *this;
92 }
93 __hostdev__ Extrema& max(const ValueT& v)
94 {
95 if (v > mMax) mMax = v;
96 return *this;
97 }
98 __hostdev__ Extrema& add(const ValueT& v)
99 {
100 this->min(v);
101 this->max(v);
102 return *this;
103 }
104 __hostdev__ Extrema& add(const ValueT& v, uint64_t) { return this->add(v); }
106 {
107 this->min(other.mMin);
108 this->max(other.mMax);
109 return *this;
110 }
111 __hostdev__ const ValueT& min() const { return mMin; }
112 __hostdev__ const ValueT& max() const { return mMax; }
113 __hostdev__ operator bool() const { return mMin <= mMax; }
114 __hostdev__ static constexpr bool hasMinMax() { return !std::is_same<bool, ValueT>::value; }
115 __hostdev__ static constexpr bool hasAverage() { return false; }
116 __hostdev__ static constexpr bool hasStdDeviation() { return false; }
117 __hostdev__ static constexpr bool hasStats() { return !std::is_same<bool, ValueT>::value; }
118 __hostdev__ static constexpr size_t size() { return 0; }
119
120 template <typename NodeT>
121 __hostdev__ void setStats(NodeT &node) const
122 {
123 node.setMin(this->min());
124 node.setMax(this->max());
125 }
126}; // Extrema<T, 0>
127
128/// @brief Template specialization of Extrema on vector value types, i.e. rank = 1
129template<typename VecT>
130class Extrema<VecT, 1>
131{
132protected:
133 using Real = typename VecT::ValueType; // this works with both nanovdb and openvdb vectors
134 struct Pair
135 {
137 VecT vector;
138
139 __hostdev__ Pair(Real s)// is only used by Extrema() default c-tor
140 : scalar(s)
141 , vector(s)
142 {
143 }
144 __hostdev__ Pair(const VecT& v)
145 : scalar(v.lengthSqr())
146 , vector(v)
147 {
148 }
149 __hostdev__ bool operator<(const Pair& rhs) const { return scalar < rhs.scalar; }
150 } mMin, mMax;
152 {
153 if (p < mMin) mMin = p;
154 if (mMax < p) mMax = p;
155 return *this;
156 }
157
158public:
159 using ValueType = VecT;
161#if defined(__CUDACC__)
162 : mMin(cuda::std::numeric_limits<Real>::max())
163 , mMax(cuda::std::numeric_limits<Real>::lowest())
164#else
165 : mMin(std::numeric_limits<Real>::max())
166 , mMax(std::numeric_limits<Real>::lowest())
167#endif
168 {
169 }
170 __hostdev__ Extrema(const VecT& v)
171 : mMin(v)
172 , mMax(v)
173 {
174 }
175 __hostdev__ Extrema(const VecT& a, const VecT& b)
176 : mMin(a)
177 , mMax(b)
178 {
179 }
180 __hostdev__ Extrema& min(const VecT& v)
181 {
182 Pair tmp(v);
183 if (tmp < mMin) mMin = tmp;
184 return *this;
185 }
186 __hostdev__ Extrema& max(const VecT& v)
187 {
188 Pair tmp(v);
189 if (mMax < tmp) mMax = tmp;
190 return *this;
191 }
192 __hostdev__ Extrema& add(const VecT& v) { return this->add(Pair(v)); }
193 __hostdev__ Extrema& add(const VecT& v, uint64_t) { return this->add(Pair(v)); }
195 {
196 if (other.mMin < mMin) mMin = other.mMin;
197 if (mMax < other.mMax) mMax = other.mMax;
198 return *this;
199 }
200 __hostdev__ const VecT& min() const { return mMin.vector; }
201 __hostdev__ const VecT& max() const { return mMax.vector; }
202 __hostdev__ operator bool() const { return !(mMax < mMin); }
203 __hostdev__ static constexpr bool hasMinMax() { return !std::is_same<bool, Real>::value; }
204 __hostdev__ static constexpr bool hasAverage() { return false; }
205 __hostdev__ static constexpr bool hasStdDeviation() { return false; }
206 __hostdev__ static constexpr bool hasStats() { return !std::is_same<bool, Real>::value; }
207 __hostdev__ static constexpr size_t size() { return 0; }
208
209 template <typename NodeT>
210 __hostdev__ void setStats(NodeT &node) const
211 {
212 node.setMin(this->min());
213 node.setMax(this->max());
214 }
215}; // Extrema<T, 1>
216
217//================================================================================================
218
219template<typename ValueT, int Rank = TensorTraits<ValueT>::Rank>
220class Stats;
221
222/// @brief This class computes statistics (minimum value, maximum
223/// value, mean, variance and standard deviation) of a population
224/// of floating-point values.
225///
226/// @details variance = Mean[ (X-Mean[X])^2 ] = Mean[X^2] - Mean[X]^2,
227/// standard deviation = sqrt(variance)
228///
229/// @note This class employs incremental computation and double precision.
230template<typename ValueT>
231class Stats<ValueT, 0> : public Extrema<ValueT, 0>
232{
233protected:
235 using RealT = double; // for accuracy the internal precission must be 64 bit floats
236 size_t mSize;
237 double mAvg, mAux;
238
239public:
240 using ValueType = ValueT;
242 : BaseT()
243 , mSize(0)
244 , mAvg(0.0)
245 , mAux(0.0)
246 {
247 }
248 __hostdev__ Stats(const ValueT& val)
249 : BaseT(val)
250 , mSize(1)
251 , mAvg(RealT(val))
252 , mAux(0.0)
253 {
254 }
255 /// @brief Add a single sample
256 __hostdev__ Stats& add(const ValueT& val)
257 {
258 BaseT::add(val);
259 mSize += 1;
260 const double delta = double(val) - mAvg;
261 mAvg += delta / double(mSize);
262 mAux += delta * (double(val) - mAvg);
263 return *this;
264 }
265 /// @brief Add @a n samples with constant value @a val.
266 __hostdev__ Stats& add(const ValueT& val, uint64_t n)
267 {
268 const double denom = 1.0 / double(mSize + n);
269 const double delta = double(val) - mAvg;
270 mAvg += denom * delta * double(n);
271 mAux += denom * delta * delta * double(mSize) * double(n);
272 BaseT::add(val);
273 mSize += n;
274 return *this;
275 }
276
277 /// Add the samples from the other Stats instance.
278 __hostdev__ Stats& add(const Stats& other)
279 {
280 if (other.mSize > 0) {
281 const double denom = 1.0 / double(mSize + other.mSize);
282 const double delta = other.mAvg - mAvg;
283 mAvg += denom * delta * double(other.mSize);
284 mAux += other.mAux + denom * delta * delta * double(mSize) * double(other.mSize);
285 BaseT::add(other);
286 mSize += other.mSize;
287 }
288 return *this;
289 }
290
291 __hostdev__ static constexpr bool hasMinMax() { return !std::is_same<bool, ValueT>::value; }
292 __hostdev__ static constexpr bool hasAverage() { return !std::is_same<bool, ValueT>::value; }
293 __hostdev__ static constexpr bool hasStdDeviation() { return !std::is_same<bool, ValueT>::value; }
294 __hostdev__ static constexpr bool hasStats() { return !std::is_same<bool, ValueT>::value; }
295
296 __hostdev__ size_t size() const { return mSize; }
297
298 //@{
299 /// Return the arithmetic mean, i.e. average, value.
300 __hostdev__ double avg() const { return mAvg; }
301 __hostdev__ double mean() const { return mAvg; }
302 //@}
303
304 //@{
305 /// @brief Return the population variance.
306 ///
307 /// @note The unbiased sample variance = population variance * num/(num-1)
308 __hostdev__ double var() const { return mSize < 2 ? 0.0 : mAux / double(mSize); }
309 __hostdev__ double variance() const { return this->var(); }
310 //@}
311
312 //@{
313 /// @brief Return the standard deviation (=Sqrt(variance)) as
314 /// defined from the (biased) population variance.
315 __hostdev__ double std() const { return sqrt(this->var()); }
316 __hostdev__ double stdDev() const { return this->std(); }
317 //@}
318
319 template <typename NodeT>
320 __hostdev__ void setStats(NodeT &node) const
321 {
322 node.setMin(this->min());
323 node.setMax(this->max());
324 node.setAvg(this->avg());
325 node.setDev(this->std());
326 }
327}; // end Stats<T, 0>
328
329/// @brief This class computes statistics (minimum value, maximum
330/// value, mean, variance and standard deviation) of a population
331/// of floating-point values.
332///
333/// @details variance = Mean[ (X-Mean[X])^2 ] = Mean[X^2] - Mean[X]^2,
334/// standard deviation = sqrt(variance)
335///
336/// @note This class employs incremental computation and double precision.
337template<typename ValueT>
338class Stats<ValueT, 1> : public Extrema<ValueT, 1>
339{
340protected:
342 using RealT = double; // for accuracy the internal precision must be 64 bit floats
343 size_t mSize;
344 double mAvg, mAux;
345
346public:
347 using ValueType = ValueT;
349 : BaseT()
350 , mSize(0)
351 , mAvg(0.0)
352 , mAux(0.0)
353 {
354 }
355 /// @brief Add a single sample
356 __hostdev__ Stats& add(const ValueT& val)
357 {
358 typename BaseT::Pair tmp(val);
359 BaseT::add(tmp);
360 mSize += 1;
361 const double delta = tmp.scalar - mAvg;
362 mAvg += delta / double(mSize);
363 mAux += delta * (tmp.scalar - mAvg);
364 return *this;
365 }
366 /// @brief Add @a n samples with constant value @a val.
367 __hostdev__ Stats& add(const ValueT& val, uint64_t n)
368 {
369 typename BaseT::Pair tmp(val);
370 const double denom = 1.0 / double(mSize + n);
371 const double delta = tmp.scalar - mAvg;
372 mAvg += denom * delta * double(n);
373 mAux += denom * delta * delta * double(mSize) * double(n);
374 BaseT::add(tmp);
375 mSize += n;
376 return *this;
377 }
378
379 /// Add the samples from the other Stats instance.
380 __hostdev__ Stats& add(const Stats& other)
381 {
382 if (other.mSize > 0) {
383 const double denom = 1.0 / double(mSize + other.mSize);
384 const double delta = other.mAvg - mAvg;
385 mAvg += denom * delta * double(other.mSize);
386 mAux += other.mAux + denom * delta * delta * double(mSize) * double(other.mSize);
387 BaseT::add(other);
388 mSize += other.mSize;
389 }
390 return *this;
391 }
392
393 __hostdev__ static constexpr bool hasMinMax() { return !std::is_same<bool, ValueT>::value; }
394 __hostdev__ static constexpr bool hasAverage() { return !std::is_same<bool, ValueT>::value; }
395 __hostdev__ static constexpr bool hasStdDeviation() { return !std::is_same<bool, ValueT>::value; }
396 __hostdev__ static constexpr bool hasStats() { return !std::is_same<bool, ValueT>::value; }
397
398 __hostdev__ size_t size() const { return mSize; }
399
400 //@{
401 /// Return the arithmetic mean, i.e. average, value.
402 __hostdev__ double avg() const { return mAvg; }
403 __hostdev__ double mean() const { return mAvg; }
404 //@}
405
406 //@{
407 /// @brief Return the population variance.
408 ///
409 /// @note The unbiased sample variance = population variance * num/(num-1)
410 __hostdev__ double var() const { return mSize < 2 ? 0.0 : mAux / double(mSize); }
411 __hostdev__ double variance() const { return this->var(); }
412 //@}
413
414 //@{
415 /// @brief Return the standard deviation (=Sqrt(variance)) as
416 /// defined from the (biased) population variance.
417 __hostdev__ double std() const { return sqrt(this->var()); }
418 __hostdev__ double stdDev() const { return this->std(); }
419 //@}
420
421 template <typename NodeT>
422 __hostdev__ void setStats(NodeT &node) const
423 {
424 node.setMin(this->min());
425 node.setMax(this->max());
426 node.setAvg(this->avg());
427 node.setDev(this->std());
428 }
429}; // end Stats<T, 1>
430
431/// @brief No-op Stats class
432template<typename ValueT>
434{
435 using ValueType = ValueT;
437 __hostdev__ NoopStats(const ValueT&) {}
438 __hostdev__ NoopStats& add(const ValueT&) { return *this; }
439 __hostdev__ NoopStats& add(const ValueT&, uint64_t) { return *this; }
440 __hostdev__ NoopStats& add(const NoopStats&) { return *this; }
441 __hostdev__ static constexpr size_t size() { return 0; }
442 __hostdev__ static constexpr bool hasMinMax() { return false; }
443 __hostdev__ static constexpr bool hasAverage() { return false; }
444 __hostdev__ static constexpr bool hasStdDeviation() { return false; }
445 __hostdev__ static constexpr bool hasStats() { return false; }
446 template <typename NodeT>
447 __hostdev__ void setStats(NodeT&) const{}
448}; // end NoopStats<T>
449
450//================================================================================================
451
452/// @brief Allows for the construction of NanoVDB grids without any dependency
453template<typename GridT, typename StatsT = Stats<typename GridT::ValueType>>
455{
456 struct NodeStats;
457 using TreeT = typename GridT::TreeType;
458 using ValueT = typename TreeT::ValueType;
459 using BuildT = typename TreeT::BuildType;
460 using Node0 = typename TreeT::Node0; // leaf
461 using Node1 = typename TreeT::Node1; // lower
462 using Node2 = typename TreeT::Node2; // upper
463 using RootT = typename TreeT::Node3; // root
464 static_assert(std::is_same<ValueT, typename StatsT::ValueType>::value, "Mismatching type");
465
466 ValueT mDelta; // skip rendering of node if: node.max < -mDelta || node.min > mDelta
467
468 void process( GridT& );// process grid and all tree nodes
469 void process( TreeT& );// process Tree, root node and child nodes
470 void process( RootT& );// process root node and child nodes
471 NodeStats process( Node0& );// process leaf node
472
473 template<typename NodeT>
474 NodeStats process( NodeT& );// process internal node and child nodes
475
476 template<typename DataT, int Rank>
477 void setStats(DataT*, const Extrema<ValueT, Rank>&);
478 template<typename DataT, int Rank>
479 void setStats(DataT*, const Stats<ValueT, Rank>&);
480 template<typename DataT>
481 void setStats(DataT*, const NoopStats<ValueT>&) {}
482
483 template<typename T, typename FlagT>
484 typename std::enable_if<!std::is_floating_point<T>::value>::type
485 setFlag(const T&, const T&, FlagT& flag) const { flag &= ~FlagT(1); } // unset 1st bit to enable rendering
486
487 template<typename T, typename FlagT>
488 typename std::enable_if<std::is_floating_point<T>::value>::type
489 setFlag(const T& min, const T& max, FlagT& flag) const;
490
491public:
492 GridStats() = default;
493
494 void operator()(GridT& grid, ValueT delta = ValueT(0));
495
496}; // GridStats
497
498template<typename GridT, typename StatsT>
499struct GridStats<GridT, StatsT>::NodeStats
500{
501 StatsT stats;
503
504 NodeStats(): stats(), bbox() {}//activeCount(0), bbox() {};
505
506 NodeStats& add(const NodeStats &other)
507 {
508 stats.add( other.stats );// no-op for NoopStats?!
509 bbox[0].minComponent(other.bbox[0]);
510 bbox[1].maxComponent(other.bbox[1]);
511 return *this;
512 }
513};// GridStats::NodeStats
514
515//================================================================================================
516
517template<typename GridT, typename StatsT>
518void GridStats<GridT, StatsT>::operator()(GridT& grid, ValueT delta)
519{
520 mDelta = delta; // delta = voxel size for level sets, else 0
521 this->process( grid );
522}
523
524//================================================================================================
525
526template<typename GridT, typename StatsT>
527template<typename DataT, int Rank>
529 setStats(DataT* data, const Extrema<ValueT, Rank>& e)
530{
531 data->setMin(e.min());
532 data->setMax(e.max());
533}
534
535template<typename GridT, typename StatsT>
536template<typename DataT, int Rank>
537inline void GridStats<GridT, StatsT>::
538 setStats(DataT* data, const Stats<ValueT, Rank>& s)
539{
540 data->setMin(s.min());
541 data->setMax(s.max());
542 data->setAvg(s.avg());
543 data->setDev(s.std());
544}
545
546//================================================================================================
547
548template<typename GridT, typename StatsT>
549template<typename T, typename FlagT>
550inline typename std::enable_if<std::is_floating_point<T>::value>::type
551GridStats<GridT, StatsT>::
552 setFlag(const T& min, const T& max, FlagT& flag) const
553{
554 if (mDelta > 0 && (min > mDelta || max < -mDelta)) {// LS: min > dx || max < -dx
555 flag |= FlagT(1u);// set 1st bit to disable rendering
556 } else {
557 flag &= ~FlagT(1u);// unset 1st bit to enable rendering
558 }
559}
560
561//================================================================================================
562
563template<typename GridT, typename StatsT>
564void GridStats<GridT, StatsT>::process( GridT &grid )
565{
566 this->process( grid.tree() );// this processes tree, root and all nodes
567
568 // set world space AABB
569 auto& data = *grid.data();
570 const auto& indexBBox = grid.tree().root().bbox();
571 if (indexBBox.empty()) {
572 data.mWorldBBox = BBox<Vec3d>();
573 data.setBBoxOn(false);
574 } else {
575 // Note that below max is offset by one since CoordBBox.max is inclusive
576 // while bbox<Vec3d>.max is exclusive. However, min is inclusive in both
577 // CoordBBox and BBox<Vec3d>. This also guarantees that a grid with a single
578 // active voxel, does not have an empty world bbox! E.g. if a grid with a
579 // unit index-to-world transformation only contains the active voxel (0,0,0)
580 // then indeBBox = (0,0,0) -> (0,0,0) and then worldBBox = (0.0, 0.0, 0.0)
581 // -> (1.0, 1.0, 1.0). This is a consequence of the different definitions
582 // of index and world bounding boxes inherited from OpenVDB!
583 grid.mWorldBBox = CoordBBox(indexBBox[0], indexBBox[1].offsetBy(1)).transform(grid.map());
584 grid.setBBoxOn(true);
585 }
586
587 // set bit flags
588 data.setMinMaxOn(StatsT::hasMinMax());
589 data.setAverageOn(StatsT::hasAverage());
590 data.setStdDeviationOn(StatsT::hasStdDeviation());
591} // GridStats::process( Grid )
592
593//================================================================================================
594
595template<typename GridT, typename StatsT>
596inline void GridStats<GridT, StatsT>::process( typename GridT::TreeType &tree )
597{
598 this->process( tree.root() );
599}
600
601//================================================================================================
602
603template<typename GridT, typename StatsT>
604void GridStats<GridT, StatsT>::process(RootT &root)
605{
606 using ChildT = Node2;
607 auto &data = *root.data();
608 if (data.mTableSize == 0) { // empty root node
609 data.mMinimum = data.mMaximum = data.mBackground;
610 data.mAverage = data.mStdDevi = 0;
611 data.mBBox = CoordBBox();
612 } else {
613 NodeStats total;
614 for (uint32_t i = 0; i < data.mTableSize; ++i) {
615 auto* tile = data.tile(i);
616 if (tile->isChild()) { // process child node
617 total.add( this->process( *data.getChild(tile) ) );
618 } else if (tile->state) { // active tile
619 const Coord ijk = tile->origin();
620 total.bbox[0].minComponent(ijk);
621 total.bbox[1].maxComponent(ijk + Coord(ChildT::DIM - 1));
622 if (StatsT::hasStats()) { // resolved at compile time
623 total.stats.add(tile->value, ChildT::NUM_VALUES);
624 }
625 }
626 }
627 this->setStats(&data, total.stats);
628 if (total.bbox.empty()) {
629 std::cerr << "\nWarning in GridStats: input tree only contained inactive root tiles!"
630 << "\nWhile not strictly an error it's rather suspicious!\n";
631 }
632 data.mBBox = total.bbox;
633 }
634} // GridStats::process( RootNode )
635
636//================================================================================================
637
638template<typename GridT, typename StatsT>
639template<typename NodeT>
640typename GridStats<GridT, StatsT>::NodeStats
641GridStats<GridT, StatsT>::process(NodeT &node)
642{
643 static_assert(is_same<NodeT,Node1>::value || is_same<NodeT,Node2>::value, "Incorrect node type");
644 using ChildT = typename NodeT::ChildNodeType;
645
646 NodeStats total;
647 auto* data = node.data();
648
649 // Serial processing of active tiles
650 if (const auto tileCount = data->mValueMask.countOn()) {
651 //total.activeCount = tileCount * ChildT::NUM_VALUES; // active tiles
652 for (auto it = data->mValueMask.beginOn(); it; ++it) {
653 if (StatsT::hasStats()) { // resolved at compile time
654 total.stats.add( data->mTable[*it].value, ChildT::NUM_VALUES );
655 }
656 const Coord ijk = node.offsetToGlobalCoord(*it);
657 total.bbox[0].minComponent(ijk);
658 total.bbox[1].maxComponent(ijk + Coord(int32_t(ChildT::DIM) - 1));
659 }
660 }
661
662 // Serial or parallel processing of child nodes
663 if (const size_t childCount = data->mChildMask.countOn()) {
664#ifndef NANOVDB_USE_TBB
665 for (auto it = data->mChildMask.beginOn(); it; ++it) {
666 total.add( this->process( *data->getChild(*it) ) );
667 }
668#else
669 std::unique_ptr<ChildT*[]> childNodes(new ChildT*[childCount]);
670 ChildT **ptr = childNodes.get();
671 for (auto it = data->mChildMask.beginOn(); it; ++it) {
672 *ptr++ = data->getChild( *it );
673 }
674 using RangeT = tbb::blocked_range<size_t>;
675 total.add( tbb::parallel_reduce(RangeT(0, childCount), NodeStats(),
676 [&](const RangeT &r, NodeStats local)->NodeStats {
677 for(size_t i=r.begin(); i!=r.end(); ++i){
678 local.add( this->process( *childNodes[i] ) );
679 }
680 return local;},
681 [](NodeStats a, const NodeStats &b)->NodeStats { return a.add( b ); }
682 ));
683#endif
684 }
685
686 data->mBBox = total.bbox;
687 if (total.bbox.empty()) {
688 data->mFlags |= uint32_t(1); // set 1st bit on to disable rendering of node
689 data->mFlags &= ~uint32_t(2); // set 2nd bit off since node does not contain active values
690 } else {
691 data->mFlags |= uint32_t(2); // set 2nd bit on since node contains active values
692 if (StatsT::hasStats()) { // resolved at compile time
693 this->setStats(data, total.stats);
694 this->setFlag(data->mMinimum, data->mMaximum, data->mFlags);
695 }
696 }
697 return total;
698} // GridStats::process( InternalNode )
699
700//================================================================================================
701
702template<typename GridT, typename StatsT>
703typename GridStats<GridT, StatsT>::NodeStats
704GridStats<GridT, StatsT>::process(Node0 &leaf)
705{
706 NodeStats local;
707 if (leaf.updateBBox()) {// optionally update active bounding box (updates data->mFlags)
708 local.bbox[0] = local.bbox[1] = leaf.mBBoxMin;
709 local.bbox[1] += Coord(leaf.mBBoxDif[0], leaf.mBBoxDif[1], leaf.mBBoxDif[2]);
710 if (StatsT::hasStats()) {// resolved at compile time
711 for (auto it = leaf.cbeginValueOn(); it; ++it) local.stats.add(*it);
712 this->setStats(&leaf, local.stats);
713 this->setFlag(leaf.getMin(), leaf.getMax(), leaf.mFlags);
714 }
715 }
716 return local;
717} // GridStats::process( LeafNode )
718
719//================================================================================================
720
721template<typename BuildT>
723{
724 using GridT = NanoGrid<BuildT>;
725 using ValueT = typename GridT::ValueType;
726 if (mode == StatsMode::Disable) {
727 return;
728 } else if (mode == StatsMode::BBox || std::is_same<bool, ValueT>::value) {
730 stats(grid);
731 } else if (mode == StatsMode::MinMax) {
733 stats(grid);
734 } else if (mode == StatsMode::All) {
736 stats(grid);
737 } else {
738 throw std::runtime_error("gridStats: Unsupported statistics mode.");
739 }
740}// gridStats
741
742//================================================================================================
743
744namespace {
745
746// returns a bitmask (of size 32^3 or 16^3) that marks all the entries
747// in a node table that intersects with the specified bounding box.
748template<typename NodeT>
749Mask<NodeT::LOG2DIM> getBBoxMask(const CoordBBox &bbox, const NodeT* node)
750{
751 Mask<NodeT::LOG2DIM> mask;// typically 32^3 or 16^3 bit mask
752 auto b = CoordBBox::createCube(node->origin(), node->dim());
753 assert( bbox.hasOverlap(b) );
754 if ( bbox.isInside(b) ) {
755 mask.setOn();//node is completely inside the bbox so early out
756 } else {
757 b.intersect(bbox);// trim bounding box
758 // transform bounding box from global to local coordinates
759 b.min() &= NodeT::DIM-1u;
760 b.min() >>= NodeT::ChildNodeType::TOTAL;
761 b.max() &= NodeT::DIM-1u;
762 b.max() >>= NodeT::ChildNodeType::TOTAL;
763 assert( !b.empty() );
764 auto it = b.begin();// iterates over all the child nodes or tiles that intersects bbox
765 for (const Coord& ijk = *it; it; ++it) {
766 mask.setOn(ijk[2] + (ijk[1] << NodeT::LOG2DIM) + (ijk[0] << 2*NodeT::LOG2DIM));
767 }
768 }
769 return mask;
770}// getBBoxMask
771
772}// end of unnamed namespace
773
774/// @brief return the extrema of all the values in a grid that
775/// intersects the specified bounding box.
776template<typename BuildT>
777Extrema<typename NanoGrid<BuildT>::ValueType>
778getExtrema(const NanoGrid<BuildT>& grid, const CoordBBox &bbox)
779{
780 using GridT = NanoGrid<BuildT>;
781 using ValueT = typename GridT::ValueType;
782 using TreeT = typename GridTree<GridT>::type;
783 using RootT = typename NodeTrait<TreeT, 3>::type;// root node
784 using Node2 = typename NodeTrait<TreeT, 2>::type;// upper internal node
785 using Node1 = typename NodeTrait<TreeT, 1>::type;// lower internal node
786 using Node0 = typename NodeTrait<TreeT, 0>::type;// leaf node
787
788 Extrema<ValueT> extrema;
789 const RootT &root = grid.tree().root();
790 const auto &bbox3 = root.bbox();
791 if (bbox.isInside(bbox3)) {// bbox3 is contained inside bbox
792 extrema.min(root.minimum());
793 extrema.max(root.maximum());
794 extrema.add(root.background());
795 } else if (bbox.hasOverlap(bbox3)) {
796 const auto *data3 = root.data();
797 for (uint32_t i=0; i<data3->mTableSize; ++i) {
798 const auto *tile = data3->tile(i);
799 CoordBBox bbox2 = CoordBBox::createCube(tile->origin(), Node2::dim());
800 if (!bbox.hasOverlap(bbox2)) continue;
801 if (tile->isChild()) {
802 const Node2 *node2 = data3->getChild(tile);
803 if (bbox.isInside(bbox2)) {
804 extrema.min(node2->minimum());
805 extrema.max(node2->maximum());
806 } else {// partial intersections at level 2
807 auto *data2 = node2->data();
808 const auto bboxMask2 = getBBoxMask(bbox, node2);
809 for (auto it2 = bboxMask2.beginOn(); it2; ++it2) {
810 if (data2->mChildMask.isOn(*it2)) {
811 const Node1* node1 = data2->getChild(*it2);
812 CoordBBox bbox1 = CoordBBox::createCube(node1->origin(), Node1::dim());
813 if (bbox.isInside(bbox1)) {
814 extrema.min(node1->minimum());
815 extrema.max(node1->maximum());
816 } else {// partial intersection at level 1
817 auto *data1 = node1->data();
818 const auto bboxMask1 = getBBoxMask(bbox, node1);
819 for (auto it1 = bboxMask1.beginOn(); it1; ++it1) {
820 if (data1->mChildMask.isOn(*it1)) {
821 const Node0* node0 = data1->getChild(*it1);
822 CoordBBox bbox0 = CoordBBox::createCube(node0->origin(), Node0::dim());
823 if (bbox.isInside(bbox0)) {
824 extrema.min(node0->minimum());
825 extrema.max(node0->maximum());
826 } else {// partial intersection at level 0
827 auto *data0 = node0->data();
828 const auto bboxMask0 = getBBoxMask(bbox, node0);
829 for (auto it0 = bboxMask0.beginOn(); it0; ++it0) {
830 extrema.add(data0->getValue(*it0));
831 }
832 }// end partial intersection at level 0
833 } else {// tile at level 1
834 extrema.add(data1->mTable[*it1].value);
835 }
836 }
837 }// end of partial intersection at level 1
838 } else {// tile at level 2
839 extrema.add(data2->mTable[*it2].value);
840 }
841 }// loop over tiles and nodes at level 2
842 }// end of partial intersection at level 1
843 } else {// tile at root level
844 extrema.add(tile->value);
845 }
846 }// loop over root table
847 } else {// bbox does not overlap the grid
848 extrema.add(root.background());
849 }
850 return extrema;
851}// getExtrema
852
853} // namespace nanovdb
854
855#endif // NANOVDB_GRIDSTATS_H_HAS_BEEN_INCLUDED
A unified wrapper for tbb::parallel_for and a naive std::thread fallback.
Implements a light-weight self-contained VDB data-structure in a single file! In other words,...
#define __hostdev__
Definition NanoVDB.h:213
Custom Range class that is compatible with the tbb::blocked_range classes.
Template specialization of Extrema on scalar value types, i.e. rank = 0.
Definition GridStats.h:62
__hostdev__ void setStats(NodeT &node) const
Definition GridStats.h:121
static __hostdev__ constexpr bool hasAverage()
Definition GridStats.h:115
static __hostdev__ constexpr size_t size()
Definition GridStats.h:118
__hostdev__ Extrema & add(const ValueT &v, uint64_t)
Definition GridStats.h:104
__hostdev__ Extrema & max(const ValueT &v)
Definition GridStats.h:93
__hostdev__ const ValueT & min() const
Definition GridStats.h:111
ValueT ValueType
Definition GridStats.h:67
static __hostdev__ constexpr bool hasMinMax()
Definition GridStats.h:114
__hostdev__ Extrema()
Definition GridStats.h:68
__hostdev__ const ValueT & max() const
Definition GridStats.h:112
__hostdev__ Extrema & add(const Extrema &other)
Definition GridStats.h:105
ValueT mMax
Definition GridStats.h:64
static __hostdev__ constexpr bool hasStdDeviation()
Definition GridStats.h:116
__hostdev__ Extrema & add(const ValueT &v)
Definition GridStats.h:98
__hostdev__ Extrema(const ValueT &a, const ValueT &b)
Definition GridStats.h:83
__hostdev__ Extrema & min(const ValueT &v)
Definition GridStats.h:88
static __hostdev__ constexpr bool hasStats()
Definition GridStats.h:117
__hostdev__ Extrema(const ValueT &v)
Definition GridStats.h:78
__hostdev__ void setStats(NodeT &node) const
Definition GridStats.h:210
static __hostdev__ constexpr bool hasAverage()
Definition GridStats.h:204
static __hostdev__ constexpr size_t size()
Definition GridStats.h:207
__hostdev__ Extrema & min(const VecT &v)
Definition GridStats.h:180
__hostdev__ Extrema & add(const Pair &p)
Definition GridStats.h:151
static __hostdev__ constexpr bool hasMinMax()
Definition GridStats.h:203
typename VecT::ValueType Real
Definition GridStats.h:133
__hostdev__ Extrema()
Definition GridStats.h:160
__hostdev__ Extrema & add(const VecT &v, uint64_t)
Definition GridStats.h:193
__hostdev__ const VecT & max() const
Definition GridStats.h:201
__hostdev__ Extrema & max(const VecT &v)
Definition GridStats.h:186
__hostdev__ Extrema(const VecT &v)
Definition GridStats.h:170
__hostdev__ Extrema & add(const Extrema &other)
Definition GridStats.h:194
__hostdev__ Extrema(const VecT &a, const VecT &b)
Definition GridStats.h:175
__hostdev__ const VecT & min() const
Definition GridStats.h:200
static __hostdev__ constexpr bool hasStdDeviation()
Definition GridStats.h:205
__hostdev__ Extrema & add(const VecT &v)
Definition GridStats.h:192
static __hostdev__ constexpr bool hasStats()
Definition GridStats.h:206
VecT ValueType
Definition GridStats.h:159
Definition GridStats.h:57
Allows for the construction of NanoVDB grids without any dependency.
Definition GridStats.h:455
void operator()(GridT &grid, ValueT delta=ValueT(0))
Definition GridStats.h:518
Highest level of the data structure. Contains a tree and a world->index transform (that currently onl...
Definition NanoVDB.h:3699
const TreeT & tree() const
Return a const reference to the tree.
Definition NanoVDB.h:3753
__hostdev__ Stats & add(const Stats &other)
Add the samples from the other Stats instance.
Definition GridStats.h:278
__hostdev__ void setStats(NodeT &node) const
Definition GridStats.h:320
static __hostdev__ constexpr bool hasAverage()
Definition GridStats.h:292
__hostdev__ double stdDev() const
Definition GridStats.h:316
__hostdev__ Stats()
Definition GridStats.h:241
ValueT ValueType
Definition GridStats.h:240
__hostdev__ double mean() const
Definition GridStats.h:301
static __hostdev__ constexpr bool hasMinMax()
Definition GridStats.h:291
size_t mSize
Definition GridStats.h:236
double RealT
Definition GridStats.h:235
__hostdev__ Stats & add(const ValueT &val)
Add a single sample.
Definition GridStats.h:256
__hostdev__ double var() const
Return the population variance.
Definition GridStats.h:308
__hostdev__ Stats & add(const ValueT &val, uint64_t n)
Add n samples with constant value val.
Definition GridStats.h:266
__hostdev__ double avg() const
Return the arithmetic mean, i.e. average, value.
Definition GridStats.h:300
__hostdev__ double std() const
Return the standard deviation (=Sqrt(variance)) as defined from the (biased) population variance.
Definition GridStats.h:315
__hostdev__ Stats(const ValueT &val)
Definition GridStats.h:248
__hostdev__ size_t size() const
Definition GridStats.h:296
static __hostdev__ constexpr bool hasStdDeviation()
Definition GridStats.h:293
double mAux
Definition GridStats.h:237
__hostdev__ double variance() const
Definition GridStats.h:309
static __hostdev__ constexpr bool hasStats()
Definition GridStats.h:294
__hostdev__ Stats & add(const Stats &other)
Add the samples from the other Stats instance.
Definition GridStats.h:380
__hostdev__ void setStats(NodeT &node) const
Definition GridStats.h:422
static __hostdev__ constexpr bool hasAverage()
Definition GridStats.h:394
__hostdev__ double stdDev() const
Definition GridStats.h:418
__hostdev__ Stats()
Definition GridStats.h:348
ValueT ValueType
Definition GridStats.h:347
__hostdev__ double mean() const
Definition GridStats.h:403
static __hostdev__ constexpr bool hasMinMax()
Definition GridStats.h:393
size_t mSize
Definition GridStats.h:343
double RealT
Definition GridStats.h:342
__hostdev__ Stats & add(const ValueT &val)
Add a single sample.
Definition GridStats.h:356
__hostdev__ double var() const
Return the population variance.
Definition GridStats.h:410
__hostdev__ Stats & add(const ValueT &val, uint64_t n)
Add n samples with constant value val.
Definition GridStats.h:367
__hostdev__ double avg() const
Return the arithmetic mean, i.e. average, value.
Definition GridStats.h:402
__hostdev__ double std() const
Return the standard deviation (=Sqrt(variance)) as defined from the (biased) population variance.
Definition GridStats.h:417
__hostdev__ size_t size() const
Definition GridStats.h:398
static __hostdev__ constexpr bool hasStdDeviation()
Definition GridStats.h:395
double mAux
Definition GridStats.h:344
__hostdev__ double variance() const
Definition GridStats.h:411
static __hostdev__ constexpr bool hasStats()
Definition GridStats.h:396
Definition GridStats.h:220
void add(double val)
Add a single sample.
Definition Stats.h:102
double min() const
Return the minimum value.
Definition Stats.h:121
double max() const
Return the maximum value.
Definition Stats.h:124
Definition NanoVDB.h:247
Grid< NanoTree< BuildT > > NanoGrid
Definition NanoVDB.h:6449
Extrema< typename NanoGrid< BuildT >::ValueType > getExtrema(const NanoGrid< BuildT > &grid, const CoordBBox &bbox)
return the extrema of all the values in a grid that intersects the specified bounding box.
Definition GridStats.h:778
BBox< Coord > CoordBBox
Definition NanoVDB.h:2535
StatsMode
Grid flags which indicate what extra information is present in the grid buffer.
Definition GridStats.h:38
void gridStats(NanoGrid< BuildT > &grid, StatsMode mode=StatsMode::Default)
Re-computes the min/max, stats and bbox information for an existing NanoVDB Grid.
Definition GridStats.h:722
Definition Coord.h:589
Definition NanoVDB.h:2295
Definition GridStats.h:135
VecT vector
Definition GridStats.h:137
__hostdev__ Pair(Real s)
Definition GridStats.h:139
__hostdev__ bool operator<(const Pair &rhs) const
Definition GridStats.h:149
Real scalar
Definition GridStats.h:136
__hostdev__ Pair(const VecT &v)
Definition GridStats.h:144
Definition GridStats.h:500
NodeStats()
Definition GridStats.h:504
StatsT stats
Definition GridStats.h:501
CoordBBox bbox
Definition GridStats.h:502
NodeStats & add(const NodeStats &other)
Definition GridStats.h:506
typename GridT::TreeType type
Definition NanoVDB.h:3976
Struct to derive node type from its level in a given grid, tree or root while preserving constness.
Definition NanoVDB.h:3404
No-op Stats class.
Definition GridStats.h:434
static __hostdev__ constexpr bool hasAverage()
Definition GridStats.h:443
__hostdev__ NoopStats & add(const ValueT &, uint64_t)
Definition GridStats.h:439
static __hostdev__ constexpr size_t size()
Definition GridStats.h:441
__hostdev__ void setStats(NodeT &) const
Definition GridStats.h:447
ValueT ValueType
Definition GridStats.h:435
static __hostdev__ constexpr bool hasMinMax()
Definition GridStats.h:442
__hostdev__ NoopStats & add(const NoopStats &)
Definition GridStats.h:440
__hostdev__ NoopStats & add(const ValueT &)
Definition GridStats.h:438
static __hostdev__ constexpr bool hasStdDeviation()
Definition GridStats.h:444
__hostdev__ NoopStats(const ValueT &)
Definition GridStats.h:437
__hostdev__ NoopStats()
Definition GridStats.h:436
static __hostdev__ constexpr bool hasStats()
Definition GridStats.h:445
static constexpr bool value
Definition NanoVDB.h:443