CondFormats/PhysicsToolsObjects/interface/Histogram.icc - master - githubmirror/cmssw

CondFormats/PhysicsToolsObjects/interface/Histogram.icc
261 строка11 KB
Cms Build
Clang-Format
24 май 2019, 00:14
24 май 2019, 00:14ecf11e6
#ifndef CondFormats_PhysicsToolsObjects_Histogram_icc
#define CondFormats_PhysicsToolsObjects_Histogram_icc

#include <cmath>
#include <cstddef>
#include <algorithm>
#include <numeric>
#include <vector>

#include "FWCore/Utilities/interface/Exception.h"

#include "CondFormats/PhysicsToolsObjects/interface/Histogram.h"

namespace PhysicsTools {
  namespace Calibration {

    namespace {  // hide when instantiated in -O0 debug mode
      template <typename T>
      inline const T &clamp(const T &min, const T &val, const T &max) {
        if (val < min)
          return min;
        if (val > max)
          return max;
        return val;
      }
    }  // namespace

    template <typename Value_t, typename Axis_t>
    Histogram<Value_t, Axis_t>::Histogram() : limits(Axis_t(), Axis_t()), total(Value_t()) {
      totalValid.store(false, std::memory_order_release);
    }

    template <typename Value_t, typename Axis_t>
    Histogram<Value_t, Axis_t>::Histogram(const Histogram &orig)
        : binULimits(orig.binULimits), binValues(orig.binValues), limits(orig.limits), total(orig.total) {
      totalValid.store(orig.totalValid.load(std::memory_order_acquire), std::memory_order_release);
    }

    template <typename Value_t, typename Axis_t>
    template <typename OValue_t, typename OAxis_t>
    Histogram<Value_t, Axis_t>::Histogram(const Histogram<OValue_t, OAxis_t> &orig)
        : binULimits(orig.binULimits.begin(), orig.binULimits.end()),
          binValues(orig.binValues.begin(), orig.binValues.end()),
          limits(orig.limits),
          total(Value_t()) {
      totalValid.store(false, std::memory_order_release);
    }

    template <typename Value_t, typename Axis_t>
    Histogram<Value_t, Axis_t>::Histogram(const std::vector<Axis_t> &binULimits)
        : binULimits(binULimits), binValues(binULimits.size() + 1), limits(Axis_t(), Axis_t()), total(Value_t()) {
      totalValid.store(true, std::memory_order_release);
      if (binULimits.size() < 2)
        throw cms::Exception("TooFewBinsError") << "Trying to generate degenerate histogram: "
                                                   "Fewer than one bin requested."
                                                << std::endl;

      limits.min = binULimits.front();
      limits.max = binULimits.back();
    }

    template <typename Value_t, typename Axis_t>
    template <typename OAxis_t>
    Histogram<Value_t, Axis_t>::Histogram(const std::vector<OAxis_t> &binULimits)
        : binULimits(binULimits.begin(), binULimits.end()),
          binValues(binULimits.size() + 1),
          limits(Axis_t(), Axis_t()),
          total(Value_t()) {
      totalValid.store(true, std::memory_order_release);
      if (binULimits.size() < 2)
        throw cms::Exception("TooFewBinsError") << "Trying to generate degenerate histogram: "
                                                   "Fewer than one bin requested."
                                                << std::endl;

      limits.min = binULimits.front();
      limits.max = binULimits.back();
    }

    template <typename Value_t, typename Axis_t>
    template <typename OAxis_t>
    Histogram<Value_t, Axis_t>::Histogram(unsigned int nBins, const PhysicsTools::Calibration::Range<OAxis_t> &range)
        : binValues(nBins + 2), limits(range), total(Value_t()) {
      totalValid.store(true, std::memory_order_release);
      if (!nBins)
        throw cms::Exception("TooFewBinsError") << "Trying to generate degenerate histogram: "
                                                   "Fewer than one bin requested."
                                                << std::endl;
    }

    template <typename Value_t, typename Axis_t>
    Histogram<Value_t, Axis_t>::Histogram(unsigned int nBins, Axis_t min, Axis_t max)
        : binValues(nBins + 2), limits(min, max), total(Value_t()) {
      totalValid.store(true, std::memory_order_release);
      if (!nBins)
        throw cms::Exception("TooFewBinsError") << "Trying to generate degenerate histogram: "
                                                   "Fewer than one bin requested."
                                                << std::endl;
    }

    template <typename Value_t, typename Axis_t>
    Histogram<Value_t, Axis_t>::~Histogram() {}

    template <typename Value_t, typename Axis_t>
    Histogram<Value_t, Axis_t> &Histogram<Value_t, Axis_t>::operator=(const Histogram &orig) {
      binULimits = orig.binULimits;
      binValues = orig.binValues;
      limits = orig.limits;
      total = orig.total;
      totalValid.store(orig.totalValid.load(std::memory_order_acquire), std::memory_order_release);
      return *this;
    }

    template <typename Value_t, typename Axis_t>
    template <typename OValue_t, typename OAxis_t>
    Histogram<Value_t, Axis_t> &Histogram<Value_t, Axis_t>::operator=(const Histogram<OValue_t, OAxis_t> &orig) {
      binULimits = std::vector<Axis_t>(orig.binULimits.begin(), orig.binULimits.end());
      binValues = std::vector<Value_t>(orig.binValues.begin(), orig.binValues.end());
      limits = orig.limits;
      total = orig.total;
      totalValid.store(orig.totalValid.load(std::memory_order_acquire), std::memory_order_release);
      return *this;
    }

    template <typename Value_t, typename Axis_t>
    void Histogram<Value_t, Axis_t>::reset() {
      std::fill(binValues.begin(), binValues.end(), Value_t());
      total = Value_t();
      totalValid.store(true, std::memory_order_release);
    }

    template <typename Value_t, typename Axis_t>
    void Histogram<Value_t, Axis_t>::setBinContent(int bin, Value_t value) {
      if (bin < 0 || (unsigned int)bin >= binValues.size())
        throw cms::Exception("RangeError") << "Histogram bin " << bin
                                           << " out of range "
                                              "[0, "
                                           << (binValues.size() - 1) << "]." << std::endl;

      binValues[bin] = value;
      totalValid.store(false, std::memory_order_release);
    }

    template <typename Value_t, typename Axis_t>
    void Histogram<Value_t, Axis_t>::fill(Axis_t x, Value_t weight) {
      int bin = findBin(x);
      binValues[bin] += weight;
      totalValid.store(false, std::memory_order_release);
    }

    template <typename Value_t, typename Axis_t>
    void Histogram<Value_t, Axis_t>::setValues(const std::vector<Value_t> &values) {
      if (values.size() != binValues.size())
        throw cms::Exception("InvalidVectorSizeError") << "Invalid vector size while setting "
                                                          "histogram values"
                                                       << std::endl;

      binValues = values;
      totalValid.store(false, std::memory_order_release);
    }

    template <typename Value_t, typename Axis_t>
    template <typename OValue_t>
    void Histogram<Value_t, Axis_t>::setValues(const std::vector<OValue_t> &values) {
      if (values.size() != binValues.size())
        throw cms::Exception("InvalidVectorSizeError") << "Invalid vector size while setting "
                                                          "histogram values"
                                                       << std::endl;

      std::copy(values.begin(), values.end(), binValues.begin());
      totalValid.store(false, std::memory_order_release);
    }

    template <typename Value_t, typename Axis_t>
    typename Histogram<Value_t, Axis_t>::Range Histogram<Value_t, Axis_t>::binRange(int bin) const {
      std::size_t size = binValues.size();
      if (bin < 1 || (unsigned int)bin > size - 2)
        throw cms::Exception("RangeError") << "Histogram bin " << bin << " out of range "
                                           << "[1, " << (size - 2) << "]." << std::endl;

      if (hasEquidistantBins()) {
        Axis_t min = (Axis_t)(bin - 1) / (size - 2);
        Axis_t max = (Axis_t)bin / (size - 2);
        min *= limits.width();
        min += limits.min;
        max *= limits.width();
        max += limits.min;
        return Range(min, max);
      } else
        return Range(binULimits[bin - 1], binULimits[bin]);
    }

    template <typename Value_t, typename Axis_t>
    int Histogram<Value_t, Axis_t>::findBin(Axis_t x) const {
      if (hasEquidistantBins()) {
        std::size_t size = binValues.size();
        x -= limits.min;
        x *= size - 2;
        x /= limits.width();
        return clamp(0, (int)(std::floor(x)) + 1, (int)size - 1);
      } else
        return std::upper_bound(binULimits.begin(), binULimits.end(), x) - binULimits.begin();
    }

    template <typename Value_t, typename Axis_t>
    Value_t Histogram<Value_t, Axis_t>::normalization() const {
      if (!totalValid.load(std::memory_order_acquire)) {
        total = std::accumulate(binValues.begin() + 1, binValues.end() - 1, Value_t());
        totalValid.store(true, std::memory_order_release);
      }

      return total;
    }

    template <typename Value_t, typename Axis_t>
    Value_t Histogram<Value_t, Axis_t>::integral(Axis_t hBound, Axis_t lBound, int mode) const {
      if (hBound < lBound)
        throw cms::Exception("InvalidIntervalError") << "Upper boundary below lower boundary in "
                                                     << "histogram integral." << std::endl;

      std::size_t size = binValues.size();
      int lBin = clamp(1, findBin(lBound), (int)size - 2);
      int hBin = clamp(1, findBin(hBound), (int)size - 2);

      double sum = Value_t();
      Range lBinRange, hBinRange;

      if (hBin > lBin)
        sum = std::accumulate(binValues.begin() + (lBin + 1), binValues.begin() + hBin, Value_t());

      if (hasEquidistantBins()) {
        Axis_t binWidth = limits.width() / (size - 2);
        lBinRange = Range((lBin - 1) * binWidth, lBin * binWidth);
        hBinRange = Range((hBin - 1) * binWidth, hBin * binWidth);
      } else {
        lBinRange = Range(binULimits[lBin - 1], binULimits[lBin]);
        hBinRange = Range(binULimits[hBin - 1], binULimits[hBin]);
      }

      switch (mode) {
        case 0:
          break;
        case 1:
          lBound = clamp(lBinRange.min, lBound, lBinRange.max);
          hBound = clamp(hBinRange.min, hBound, hBinRange.max);
          sum += binValues[lBin] * (lBinRange.max - lBound) / lBinRange.width();
          sum += binValues[hBin] * (hBound - hBinRange.min) / hBinRange.width();
          break;
        default:
          throw cms::Exception("InvalidMode") << "Invalid mode " << mode
                                              << " in "
                                                 "Histogram::integral()"
                                              << std::endl;
      }

      return sum;
    }

  }  // namespace Calibration
}  // namespace PhysicsTools

#endif  // CondFormats_PhysicsToolsObjects_Histogram_icc