27 #ifndef TOMOGRAPHER_MHRW_BIN_ERR_H 28 #define TOMOGRAPHER_MHRW_BIN_ERR_H 36 #include <boost/math/constants/constants.hpp> 50 namespace tomo_internal {
54 template<
int NumLevels,
bool calculate>
55 struct helper_samples_size {
57 value = Eigen::Dynamic
60 template<
int NumLevels>
61 struct helper_samples_size<NumLevels,true> {
63 value = (1 << NumLevels)
107 template<
typename ValueType_,
int NumTrackValues_ = Eigen::Dynamic,
int NumLevels_ = Eigen::Dynamic,
108 bool StoreBinSums_ =
true,
typename CountIntType_ =
int>
121 static constexpr
int NumTrackValuesCTime = NumTrackValues_;
125 static constexpr
int NumLevelsCTime = NumLevels_;
129 static constexpr
int NumLevelsPlusOneCTime = (NumLevelsCTime == Eigen::Dynamic
131 : (NumLevelsCTime + 1));
154 static constexpr
int SamplesSizeCTime =
162 static constexpr
bool StoreBinSums = StoreBinSums_;
192 UNKNOWN_CONVERGENCE = 0,
257 template<
typename Params,
typename LoggerType_>
273 static constexpr
int NumTrackValuesCTime = Params::NumTrackValuesCTime;
277 static constexpr
int NumLevelsCTime = Params::NumLevelsCTime;
281 static constexpr
int NumLevelsPlusOneCTime = Params::NumLevelsPlusOneCTime;
285 static constexpr
int SamplesSizeCTime = Params::SamplesSizeCTime;
289 static constexpr
bool StoreBinSums = Params::StoreBinSums;
340 UNKNOWN_CONVERGENCE = Params::UNKNOWN_CONVERGENCE,
342 CONVERGED = Params::CONVERGED,
344 NOT_CONVERGED = Params::NOT_CONVERGED
360 SamplesArray samples;
368 CountIntType n_samples;
374 CountIntType n_flushes;
395 BinSumSqArray bin_sumsq;
415 : numTrackValues(num_track_values_),
416 numLevels(num_levels_),
417 samplesSize(1 << numLevels()),
418 samples(numTrackValues(), samplesSize()),
420 bin_sum(BinSumArray::Zero(numTrackValues())),
421 bin_sumsq(BinSumSqArray::Zero(numTrackValues(), numLevels()+1)),
426 tomographer_assert( (1<<numLevels()) == samplesSize() );
440 helper_reset_bin_sum();
441 bin_sumsq = BinSumSqArray::Zero(numTrackValues(), numLevels()+1);
442 logger.longdebug(
"BinningAnalysis::reset()",
"ready to go.");
455 template<
typename Derived>
458 const int ninbin = n_samples % samplesSize();
462 tomographer_assert(vals.rows() == numTrackValues());
463 tomographer_assert(vals.cols() == 1);
466 samples.col(ninbin) = vals;
469 helper_update_bin_sum(samples.col(ninbin));
472 if ( ninbin == samplesSize() - 1 ) {
476 logger.longdebug(
"BinningAnalysis::processNewValues()", [&](
std::ostream & str) {
477 str <<
"n_samples is now " << n_samples <<
"; flushing bins. samplesSize() = " << samplesSize();
483 for (
int level = 0; level <= numLevels(); ++level) {
485 const int binnedsize = 1 << (numLevels()-level);
487 logger.longdebug(
"BinningAnalysis::processNewValues()", [&](
std::ostream & str) {
488 str <<
"Processing binning level = " << level <<
": binnedsize="<<binnedsize
489 <<
"; n_flushes=" << n_flushes <<
"\n";
490 str <<
"\tbinned samples = \n" << samples.block(0,0,numTrackValues(),binnedsize);
493 for (
int ksample = 0; ksample < binnedsize; ++ksample) {
494 bin_sumsq.
col(level) += samples.col(ksample).cwiseProduct(samples.col(ksample));
495 if (ksample % 2 == 0 && binnedsize > 1) {
496 samples.col(ksample/2) = boost::math::constants::half<ValueType>() *
497 (samples.col(ksample) + samples.col(ksample+1));
503 logger.longdebug(
"BinningAnalysis::processNewValues()", [&](
std::ostream & str) {
504 str <<
"Flushing #" << n_flushes <<
" done. bin_sum is = \n" << bin_sum <<
"\n" 505 <<
"\tbin_sumsq is = \n" << bin_sumsq <<
"\n";
518 template<
typename CalcValType, TOMOGRAPHER_ENABLED_IF_TMPL(NumTrackValuesCTime == 1)>
530 TOMOGRAPHER_ENABLED_IF(StoreBinSums)
531 inline void helper_reset_bin_sum()
533 bin_sum.value = BinSumArray::Zero(numTrackValues());
535 TOMOGRAPHER_ENABLED_IF(!StoreBinSums)
536 inline void helper_reset_bin_sum() { }
538 TOMOGRAPHER_ENABLED_IF(StoreBinSums)
542 bin_sum.value += new_samples;
544 template<
typename Derived,
545 TOMOGRAPHER_ENABLED_IF_TMPL(!StoreBinSums)>
563 TOMOGRAPHER_ENABLED_IF(StoreBinSums)
564 inline auto getBinMeans() const
565 #ifndef TOMOGRAPHER_PARSED_BY_DOXYGEN 566 -> decltype(BinSumArray() / ValueType(n_samples))
569 return bin_sum.value / ValueType(n_samples);
579 #ifndef TOMOGRAPHER_PARSED_BY_DOXYGEN 581 bin_sumsq.cwiseQuotient(n_flushes * Tools::replicated<NumTrackValuesCTime,1>(
583 .transpose().reverse(),
590 return bin_sumsq.cwiseQuotient(n_flushes * Tools::replicated<NumTrackValuesCTime,1>(
592 .transpose().reverse(),
603 TOMOGRAPHER_ENABLED_IF(StoreBinSums)
604 inline const BinSumArray & getBinSum()
const {
return bin_sum.value; }
627 template<
typename Derived>
630 tomographer_assert(means.rows() == numTrackValues());
631 tomographer_assert(means.cols() == 1);
632 const int n_levels_plus_one = numLevels()+1;
633 const int n_track_values = numTrackValues();
640 getBinSqmeans() - Tools::replicated<1,NumLevelsPlusOneCTime>(
645 ).cwiseMax(0).cwiseQuotient(
647 Tools::replicated<NumTrackValuesCTime,1>(
649 .transpose().reverse(),
654 numTrackValues(), n_levels_plus_one,
672 template<
typename Derived>
674 tomographer_assert(means.rows() == numTrackValues());
675 tomographer_assert(means.cols() == 1);
677 bin_sumsq.col(numLevels()) / ValueType(n_flushes) - means.
cwiseProduct(means).template cast<ValueType>()
678 ).cwiseMax(0).cwiseSqrt() /
std::sqrt(ValueType(n_flushes-1));
692 TOMOGRAPHER_ENABLED_IF(StoreBinSums)
693 inline BinSumSqArray calcErrorLevels()
const {
694 BinSumArray means = getBinMeans();
695 return calcErrorLevels(means);
710 TOMOGRAPHER_ENABLED_IF(StoreBinSums)
711 inline BinSumArray calcErrorLastLevel()
const {
712 BinSumArray means = getBinMeans();
713 return calcErrorLastLevel(means);
730 Eigen::ArrayXi converged_status(numTrackValues());
732 tomographer_assert(error_levels.rows() == numTrackValues());
733 tomographer_assert(error_levels.cols() == numLevels() + 1);
735 logger.longdebug(
"BinningAnalysis::determineErrorConvergence", [&](
std::ostream & str) {
736 str <<
"error_levels = \n" << error_levels <<
"\n";
743 if (numLevels() < range-1) {
745 converged_status = Eigen::ArrayXi::Constant(numTrackValues(), UNKNOWN_CONVERGENCE);
749 converged_status = Eigen::ArrayXi::Constant(numTrackValues(), CONVERGED);
751 const auto & errors = error_levels.col(numLevels());
753 for (
int level = numLevels()+1 - range; level < numLevels(); ++level) {
755 const auto & errors_thislevel = error_levels.col(level);
757 logger.longdebug(
"BinningAnalysis::determineErrorConvergence", [&](
std::ostream & str) {
758 str <<
"About to study level " << level <<
": at this point, converged_status = \n" 759 << converged_status <<
"\nand errors_thislevel = \n" << errors_thislevel;
762 for (
int val_it = 0; val_it < numTrackValues(); ++val_it) {
763 if (errors_thislevel(val_it) >= errors(val_it) &&
764 converged_status(val_it) != NOT_CONVERGED) {
765 converged_status(val_it) = CONVERGED;
766 }
else if (errors_thislevel(val_it) < 0.824 * errors(val_it)) {
767 converged_status(val_it) = NOT_CONVERGED;
768 }
else if ((errors_thislevel(val_it) < 0.9 * errors(val_it)) &&
769 converged_status(val_it) != NOT_CONVERGED) {
770 converged_status(val_it) = UNKNOWN_CONVERGENCE;
778 logger.longdebug(
"BinningAnalysis::determineErrorConvergence", [&](
std::ostream & str) {
779 str <<
"Done. converged_status [0=UNNOWN,1=CONVERGED,2=NOT CONVERGED] = \n" 783 return converged_status;
void processNewValue(const CalcValType val)
Process a new value (if we're tracking a single function only)
Params::BinSumSqArray BinSumSqArray
Type used to store the sum of squares of values at each binning level. See BinningAnalysisParams::Bin...
Base namespace for the Tomographer project.
Group template parameters for BinningAnalysis.
BinSumSqArray calcErrorLevels(const Eigen::ArrayBase< Derived > &means) const
Calculate the error bars of samples at different binning levels.
Eigen::Array< ValueType, NumTrackValuesCTime, 1 > BinSumArray
Type used to store the sum of values.
Params::CountIntType CountIntType
Type used to count the number of samples. See BinningAnalysisParams::CountIntType.
Eigen::ArrayXi determineErrorConvergence(const Eigen::Ref< const BinSumSqArray > &error_levels) const
Attempt to determine if the error bars have converged.
void processNewValues(const Eigen::DenseBase< Derived > &vals)
Process new raw samples.
const Tools::StaticOrDynamic< int,(NumLevelsCTime==Eigen::Dynamic), NumLevelsCTime > numLevels
The number of levels in the binning analysis.
Params::ValueType ValueType
Type of the value(s) for which we are calculating error bars. See BinningAnalysisParams::ValueType.
const CwiseBinaryOp< internal::scalar_product_op< typename Derived::Scalar, typename OtherDerived::Scalar >, const Derived, const OtherDerived > cwiseProduct(const Eigen::ArrayBase< OtherDerived > &other) const
BinSumArray calcErrorLastLevel(const Eigen::ArrayBase< Derived > &means) const
Calculate the error bar of samples (from the last binning level).
ValueType_ ValueType
Type of the value which we are calculating error bars for. Also the type of the error bars themselves...
void reset()
Reset this object and start again as if freshly constructed.
Basic utilities for dealing with Eigen matrices and other types.
The error bars appear to have converged.
Params::SamplesArray SamplesArray
Type of the internal buffer for the raw samples. See BinningAnalysisParams::SamplesArray.
Params::BinSumArray BinSumArray
Type used to store the sum of values. See BinningAnalysisParams::BinSumArray.
Some C++ utilities, with a tad of C++11 tricks.
Simple binning analysis for determining error bars.
LoggerType_ LoggerType
Type of the logger we will be logging debugging messages to.
BinningAnalysis(int num_track_values_, int num_levels_, LoggerType &logger_)
Constructor.
Managing the need for specific overrides to operator new() for some types (especially Eigen types) ...
const BinSumSqArray & getBinSumsq() const
Get the raw sums of the squared values observed, at each binning level.
Eigen::Array< ValueType, NumTrackValuesCTime, NumLevelsPlusOneCTime > BinSumSqArray
Type used to store the sum of squares of values at each binning level.
VarValueDecoder< T >::RetType value(const Var &var)
Access the value of the given variable, as a C++ type.
const Tools::StaticOrDynamic< CountIntType,(SamplesSizeCTime==Eigen::Dynamic), SamplesSizeCTime > samplesSize
The size of our samples buffer. (avoid using, might change in the future)
Eigen::Array< ValueType, NumTrackValuesCTime, SamplesSizeCTime > SamplesArray
Type of the internal buffer for the raw samples.
CountIntType_ CountIntType
Type used to count the number of samples. Usually int is fine, except if you are really taking a LOT ...
CountIntType getNumFlushes() const
Return the number of times the collected samples were flushed.
const Tools::StaticOrDynamic< int,(NumTrackValuesCTime==Eigen::Dynamic), NumTrackValuesCTime > numTrackValues
The number of functions being tracked/analyzed.
Utilities for logging messages.
auto getBinSqmeans() const
Get the raw average of the squared values observed, for each binning level.