9 #ifndef OPENVDB_POINTS_STATISTICS_IMPL_HAS_BEEN_INCLUDED 10 #define OPENVDB_POINTS_STATISTICS_IMPL_HAS_BEEN_INCLUDED 19 namespace statistics_internal
24 template <
typename ValueT>
27 using ExtentT = std::pair<ValueT, ValueT>;
28 ScalarMinMax(
const ValueT& init) : mMinMax(init, init) {}
29 ScalarMinMax(
const ExtentT& init) : mMinMax(init) {}
30 inline void operator()(
const ValueT& b)
32 mMinMax.first =
std::min(mMinMax.first, b);
33 mMinMax.second =
std::max(mMinMax.second, b);
35 inline void operator()(
const ExtentT& b)
37 mMinMax.first =
std::min(mMinMax.first, b.first);
38 mMinMax.second =
std::max(mMinMax.second, b.second);
40 inline const ExtentT&
get()
const {
return mMinMax; }
47 template <
typename VecT,
bool MagResult = true>
48 struct MagnitudeExtent
49 :
public ScalarMinMax<typename ValueTraits<VecT>::ElementType>
52 using ExtentT =
typename ScalarMinMax<ElementT>::ExtentT;
53 using BaseT = ScalarMinMax<ElementT>;
54 MagnitudeExtent(
const VecT& init) : BaseT(init.lengthSqr()) {}
55 MagnitudeExtent(
const ExtentT& init) : BaseT(init) {}
56 inline void operator()(
const VecT& b) { this->BaseT::operator()(b.lengthSqr()); }
57 inline void operator()(
const ExtentT& b) { this->BaseT::operator()(b); }
62 template <
typename VecT>
63 struct MagnitudeExtent<VecT, false>
66 using ExtentT = std::pair<VecT, VecT>;
67 MagnitudeExtent(
const VecT& init)
68 : mLengths(), mMinMax(init, init) {
69 mLengths.first = init.lengthSqr();
70 mLengths.second = mLengths.first;
72 MagnitudeExtent(
const ExtentT& init)
73 : mLengths(), mMinMax(init) {
74 mLengths.first = init.first.lengthSqr();
75 mLengths.second = init.second.lengthSqr();
77 inline const ExtentT&
get()
const {
return mMinMax; }
78 inline void operator()(
const VecT& b)
80 const ElementT l = b.lengthSqr();
81 if (l < mLengths.first) {
85 else if (l > mLengths.second) {
90 inline void operator()(
const ExtentT& b)
92 ElementT l = b.first.lengthSqr();
93 if (l < mLengths.first) {
95 mMinMax.first = b.first;
97 l = b.second.lengthSqr();
98 if (l > mLengths.second) {
100 mMinMax.second = b.second;
104 std::pair<ElementT, ElementT> mLengths;
111 template <
typename VecT>
112 struct ComponentExtent
114 using ExtentT = std::pair<VecT, VecT>;
115 ComponentExtent(
const VecT& init) : mMinMax(init, init) {}
116 ComponentExtent(
const ExtentT& init) : mMinMax(init) {}
117 inline const ExtentT&
get()
const {
return mMinMax; }
118 inline void operator()(
const VecT& b)
123 inline void operator()(
const ExtentT& b)
132 template <
typename ValueT,
136 typename PointDataTreeT>
137 bool evalExtents(
const PointDataTreeT& points,
138 const std::string& attribute,
139 typename ExtentOp::ExtentT& ext,
140 const FilterT& filter,
141 typename PointDataTreeT::template ValueConverter
142 <
typename ExtentOp::ExtentT::first_type>::Type*
const minTree =
nullptr,
143 typename PointDataTreeT::template ValueConverter
144 <
typename ExtentOp::ExtentT::second_type>::Type*
const maxTree =
nullptr)
146 static_assert(std::is_base_of<TreeBase, PointDataTreeT>::value,
147 "PointDataTreeT in instantiation of evalExtents is not an openvdb Tree type");
150 typename ExtentOp::ExtentT ext;
154 tree::LeafManager<const PointDataTreeT> manager(points);
155 if (manager.leafCount() == 0)
return false;
156 const size_t idx = manager.leaf(0).attributeSet().find(attribute);
157 if (idx == AttributeSet::INVALID_POS)
return false;
160 std::vector<std::unique_ptr<typename ExtentOp::ExtentT>> values;
161 if (minTree || maxTree) values.resize(manager.leafCount());
163 const ResultType result = tbb::parallel_reduce(manager.leafRange(),
165 [idx, &filter, &values]
166 (
const auto& range, ResultType in) -> ResultType
168 for (
auto leaf = range.begin(); leaf; ++leaf) {
169 AttributeHandle<ValueT, CodecT> handle(leaf->constAttributeArray(idx));
170 if (handle.size() == 0)
continue;
172 const size_t size = handle.isUniform() ? 1 : handle.size();
173 ExtentOp
op(handle.get(0));
174 for (
size_t i = 1; i < size; ++i) {
178 if (!values.empty()) {
179 values[leaf.pos()].reset(
new typename ExtentOp::ExtentT(
op.get()));
181 if (in.data)
op(in.ext);
186 auto iter = leaf->beginIndexOn(filter);
188 ExtentOp
op(handle.get(*iter));
190 for (; iter; ++iter)
op(handle.get(*iter));
191 if (!values.empty()) {
192 values[leaf.pos()].reset(
new typename ExtentOp::ExtentT(
op.get()));
194 if (in.data)
op(in.ext);
202 [](
const ResultType& a,
const ResultType& b) -> ResultType {
203 if (!b.data)
return a;
204 if (!a.data)
return b;
205 ExtentOp
op(a.ext);
op(b.ext);
216 if (minTree || maxTree) {
217 manager.foreach([minTree, maxTree, &values]
218 (
const auto& leaf,
const size_t idx) {
219 const auto& v = values[idx];
220 if (v ==
nullptr)
return;
221 const Coord& origin = leaf.origin();
222 if (minTree) minTree->addTile(1, origin, v->first,
true);
223 if (maxTree) maxTree->addTile(1, origin, v->second,
true);
227 if (result.data) ext = result.ext;
231 template <
typename ValueT,
234 typename PointDataTreeT,
235 typename std::enable_if<ValueTraits<ValueT>::IsVec,
int>::type = 0>
236 bool evalExtents(
const PointDataTreeT& points,
237 const std::string& attribute,
240 const FilterT& filter,
241 typename PointDataTreeT::template ValueConverter<ValueT>::Type* minTree,
242 typename PointDataTreeT::template ValueConverter<ValueT>::Type* maxTree)
244 typename ComponentExtent<ValueT>::ExtentT ext;
245 const bool s = evalExtents<ValueT, CodecT, FilterT,
246 ComponentExtent<ValueT>, PointDataTreeT>
247 (points, attribute, ext, filter, minTree, maxTree);
248 if (s) min = ext.first, max = ext.second;
252 template <
typename ValueT,
255 typename PointDataTreeT,
256 typename std::enable_if<!ValueTraits<ValueT>::IsVec,
int>::type = 0>
257 bool evalExtents(
const PointDataTreeT& points,
258 const std::string& attribute,
261 const FilterT& filter,
262 typename PointDataTreeT::template ValueConverter<ValueT>::Type* minTree,
263 typename PointDataTreeT::template ValueConverter<ValueT>::Type* maxTree)
265 typename ScalarMinMax<ValueT>::ExtentT ext;
266 const bool s = evalExtents<ValueT, CodecT, FilterT,
267 ScalarMinMax<ValueT>, PointDataTreeT>
268 (points, attribute, ext, filter, minTree, maxTree);
269 if (s) min = ext.first, max = ext.second;
277 template <
typename ValueT,
280 typename PointDataTreeT>
282 const std::string& attribute,
285 const FilterT& filter,
286 typename PointDataTreeT::template ValueConverter<ValueT>::Type* minTree,
287 typename PointDataTreeT::template ValueConverter<ValueT>::Type* maxTree)
289 return statistics_internal::evalExtents<ValueT, CodecT, FilterT, PointDataTreeT>
290 (points, attribute,
min,
max, filter, minTree, maxTree);
293 template <
typename ValueT,
296 typename PointDataTreeT,
297 typename ResultTreeT>
299 const std::string& attribute,
301 const FilterT& filter,
302 typename PointDataTreeT::template ValueConverter<ResultTreeT>::Type* averageTree)
308 Sample(
const ResultT& _avg,
size_t _size) : avg(_avg), size(_size) {}
310 void add(
const ResultT& val)
313 const ResultT delta = val - avg;
314 avg = avg + (delta /
static_cast<double>(size));
317 void add(
const Sample& other)
320 const double denom = 1.0 /
static_cast<double>(size + other.size);
321 const ResultT delta = other.avg - avg;
322 avg = avg + (denom * delta *
static_cast<double>(other.size));
326 ResultT avg;
size_t size;
329 static_assert(std::is_base_of<TreeBase, PointDataTreeT>::value,
330 "PointDataTreeT in instantiation of evalAverage is not an openvdb Tree type");
331 static_assert(std::is_constructible<ResultT, ValueT>::value,
332 "Target value in points::evalAverage is not constructible from the source value type.");
335 if (manager.
leafCount() == 0)
return false;
336 const size_t idx = manager.
leaf(0).attributeSet().find(attribute);
337 if (idx == AttributeSet::INVALID_POS)
return false;
339 std::vector<std::unique_ptr<Sample>> values;
342 [idx, &filter, &values] (
const auto& range) {
343 for (
auto leaf = range.begin(); leaf; ++leaf) {
345 size_t size = handle.
size();
346 if (size == 0)
continue;
348 std::unique_ptr<Sample> S(
new Sample(ResultT(handle.get(0)), 1));
349 if (handle.isUniform()) {
350 S->avg = S->avg /
static_cast<double>(size);
354 for (
size_t i = 1; i < size; ++i) {
356 S->add(ResultT(handle.get(
Index(i))));
359 values[leaf.pos()] = std::move(S);
362 auto iter = leaf->beginIndexOn(filter);
364 std::unique_ptr<Sample> S(
new Sample(ResultT(handle.get(*iter)), 1));
366 for (; iter; ++iter, ++size) {
367 S->add(ResultT(handle.get(*iter)));
369 values[leaf.pos()] = std::move(S);
374 auto iter = values.cbegin();
375 while (iter != values.cend() && !(*iter)) ++iter;
376 if (iter == values.cend())
return false;
382 for (; iter != values.cend(); ++iter) {
383 if (*iter) S.add(**iter);
392 manager.
foreach([averageTree, &values]
393 (
const auto& leaf,
const size_t idx) {
394 const auto& S = values[idx];
395 if (S ==
nullptr)
return;
396 const Coord& origin = leaf.origin();
397 averageTree->addTile(1, origin, S->avg,
true);
404 template <
typename ValueT,
407 typename PointDataTreeT,
408 typename ResultTreeT>
410 const std::string& attribute,
412 const FilterT& filter,
413 typename PointDataTreeT::template ValueConverter<ResultTreeT>::Type* totalTree)
418 static_assert(std::is_base_of<TreeBase, PointDataTreeT>::value,
419 "PointDataTreeT in instantiation of accumulate is not an openvdb Tree type");
420 static_assert(std::is_constructible<ResultT, ValueT>::value,
421 "Target value in points::accumulate is not constructible from the source value type.");
424 if (manager.
leafCount() == 0)
return false;
425 const size_t idx = manager.
leaf(0).attributeSet().find(attribute);
426 if (idx == AttributeSet::INVALID_POS)
return false;
428 std::vector<std::unique_ptr<ResultT>> values;
431 [idx, &filter, &values](
const auto& range) {
432 for (
auto leaf = range.begin(); leaf; ++leaf) {
434 if (handle.size() == 0)
continue;
436 const size_t size = handle.
isUniform() ? 1 : handle.size();
437 auto total = ResultT(handle.get(0));
438 for (
size_t i = 1; i < size; ++i) {
440 total += ResultT(handle.get(
Index(i)));
442 values[leaf.pos()].reset(
new ResultT(total));
445 auto iter = leaf->beginIndexOn(filter);
447 auto total = ResultT(handle.get(*iter));
449 for (; iter; ++iter) total += ResultT(handle.get(*iter));
450 values[leaf.pos()].reset(
new ResultT(total));
455 auto iter = values.cbegin();
456 while (iter != values.cend() && !(*iter)) ++iter;
457 if (iter == values.cend())
return false;
459 total = **iter; ++iter;
461 if (std::is_integral<ElementT>::value) {
462 using RangeT = tbb::blocked_range<const std::unique_ptr<ResultT>*>;
464 total = tbb::parallel_reduce(RangeT(&(*iter), (&values.back())+1, 32), total,
465 [](
const RangeT& range, ResultT p) -> ResultT {
466 for (
const auto& r : range)
if (r) p += *r;
468 }, std::plus<ResultT>());
471 for (; iter != values.cend(); ++iter) {
472 if (*iter) total += (**iter);
481 manager.
foreach([totalTree, &values]
482 (
const auto& leaf,
const size_t idx) {
483 const auto& v = values[idx];
484 if (v ==
nullptr)
return;
485 const Coord& origin = leaf.origin();
486 totalTree->addTile(1, origin, *v,
true);
493 template <
typename ValueT,
496 typename PointDataTreeT>
497 std::pair<ValueT, ValueT>
499 const std::string& attribute,
500 const FilterT& filter)
502 std::pair<ValueT, ValueT> results {
503 zeroVal<ValueT>(), zeroVal<ValueT>()
505 evalMinMax<ValueT, CodecT, FilterT, PointDataTreeT>
506 (points, attribute, results.first, results.second, filter);
510 template <
typename ValueT,
513 typename PointDataTreeT>
516 const std::string& attribute,
517 const FilterT& filter)
520 ConvertedT result = zeroVal<ConvertedT>();
521 evalAverage<ValueT, CodecT, FilterT, PointDataTreeT>(points, attribute, result, filter);
525 template <
typename ValueT,
528 typename PointDataTreeT>
531 const std::string& attribute,
532 const FilterT& filter)
535 PromotedT result = zeroVal<PromotedT>();
536 accumulate<ValueT, CodecT, FilterT, PointDataTreeT>(points, attribute, result, filter);
544 #endif // OPENVDB_POINTS_STATISTICS_IMPL_HAS_BEEN_INCLUDED std::pair< ValueT, ValueT > evalMinMax(const PointDataTreeT &points, const std::string &attribute, const FilterT &filter=NullFilter())
Evaluates the minimum and maximum values of a point attribute.
Definition: PointStatisticsImpl.h:498
bool isUniform() const
Definition: AttributeArray.h:2087
ConvertElementType< ValueT, double >::Type evalAverage(const PointDataTreeT &points, const std::string &attribute, const FilterT &filter=NullFilter())
Evaluates the average value of a point attribute.
Definition: PointStatisticsImpl.h:515
Definition: AttributeArray.h:763
Index32 Index
Definition: Types.h:54
Definition: IndexIterator.h:44
size_t leafCount() const
Return the number of leaf nodes.
Definition: LeafManager.h:288
#define OPENVDB_ASSERT(X)
Definition: Assert.h:41
LeafRange leafRange(size_t grainsize=1) const
Return a TBB-compatible LeafRange.
Definition: LeafManager.h:346
LeafType & leaf(size_t leafIdx) const
Return a pointer to the leaf node at index leafIdx in the array.
Definition: LeafManager.h:319
Definition: Exceptions.h:13
Index size() const
Definition: AttributeArray.h:786
void foreach(const LeafOp &op, bool threaded=true, size_t grainSize=1)
Threaded method that applies a user-supplied functor to each leaf node in the LeafManager.
Definition: LeafManager.h:484
SubT Type
Definition: Types.h:320
This class manages a linear array of pointers to a given tree's leaf nodes, as well as optional auxil...
Definition: LeafManager.h:85
Vec2< T > maxComponent(const Vec2< T > &v1, const Vec2< T > &v2)
Return component-wise maximum of the two vectors.
Definition: Vec2.h:513
Vec2< T > minComponent(const Vec2< T > &v1, const Vec2< T > &v2)
Return component-wise minimum of the two vectors.
Definition: Vec2.h:504
PromoteType< ValueT >::Highest accumulate(const PointDataTreeT &points, const std::string &attribute, const FilterT &filter=NullFilter())
Evaluates the total value of a point attribute.
Definition: PointStatisticsImpl.h:530
#define OPENVDB_VERSION_NAME
The version namespace name for this library version.
Definition: version.h.in:121
typename T::ValueType ElementType
Definition: Types.h:302
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.h.in:218