[hist] Conversion functions from RHist to TH3*

hist/histv7util/CMakeLists.txt

@@ -3,10 +3,12 @@ ROOT_STANDARD_LIBRARY_PACKAGE(ROOTHistUtil
     ROOT/Hist/ConversionUtils.hxx
     ROOT/Hist/ConvertToTH1.hxx
     ROOT/Hist/ConvertToTH2.hxx
+    ROOT/Hist/ConvertToTH3.hxx
   SOURCES
     ConversionUtils.cxx
     ConvertToTH1.cxx
     ConvertToTH2.cxx
+    ConvertToTH3.cxx
   DEPENDENCIES
     Core
     Hist

hist/histv7util/inc/ROOT/Hist/ConvertToTH3.hxx

@@ -0,0 +1,120 @@
+/// \file
+/// \warning This is part of the %ROOT 7 prototype! It will change without notice. It might trigger earthquakes.
+/// Feedback is welcome!
+
+#ifndef ROOT_Hist_ConvertToTH3
+#define ROOT_Hist_ConvertToTH3
+
+#include <ROOT/RBinWithError.hxx>
+#include <ROOT/RHist.hxx>
+#include <ROOT/RHistEngine.hxx>
+
+class TH3C;
+class TH3S;
+class TH3I;
+class TH3L;
+class TH3F;
+class TH3D;
+
+#include <memory>
+
+namespace ROOT {
+namespace Experimental {
+namespace Hist {
+
+/// Convert a two-dimensional histogram to TH3C.
+///
+/// \copydetails ConvertToTH3I(const RHistEngine<int> &engine)
+std::unique_ptr<TH3C> ConvertToTH3C(const RHistEngine<char> &engine);
+
+/// Convert a two-dimensional histogram to TH3S.
+///
+/// \copydetails ConvertToTH3I(const RHistEngine<int> &engine)
+std::unique_ptr<TH3S> ConvertToTH3S(const RHistEngine<short> &engine);
+
+/// Convert a two-dimensional histogram to TH3I.
+///
+/// As RHistEngine does not have global statistics, the number of entries and the total sum of weights will be unset.
+///
+/// Throws an exception if the histogram has more than one dimension.
+///
+/// \param[in] engine the RHistEngine to convert
+/// \return the converted TH3
+std::unique_ptr<TH3I> ConvertToTH3I(const RHistEngine<int> &engine);
+
+/// Convert a two-dimensional histogram to TH3L.
+///
+/// \copydetails ConvertToTH3I(const RHistEngine<int> &engine)
+std::unique_ptr<TH3L> ConvertToTH3L(const RHistEngine<long> &engine);
+
+/// Convert a two-dimensional histogram to TH3L.
+///
+/// \copydetails ConvertToTH3I(const RHistEngine<int> &engine)
+std::unique_ptr<TH3L> ConvertToTH3L(const RHistEngine<long long> &engine);
+
+/// Convert a two-dimensional histogram to TH3F.
+///
+/// \copydetails ConvertToTH3I(const RHistEngine<int> &engine)
+std::unique_ptr<TH3F> ConvertToTH3F(const RHistEngine<float> &engine);
+
+/// Convert a two-dimensional histogram to TH3D.
+///
+/// \copydetails ConvertToTH3I(const RHistEngine<int> &engine)
+std::unique_ptr<TH3D> ConvertToTH3D(const RHistEngine<double> &engine);
+
+/// Convert a two-dimensional histogram to TH3D.
+///
+/// \copydetails ConvertToTH3I(const RHistEngine<int> &engine)
+std::unique_ptr<TH3D> ConvertToTH3D(const RHistEngine<RBinWithError> &engine);
+
+/// Convert a two-dimensional histogram to TH3C.
+///
+/// \copydetails ConvertToTH3I(const RHist<int> &hist)
+std::unique_ptr<TH3C> ConvertToTH3C(const RHist<char> &hist);
+
+/// Convert a two-dimensional histogram to TH3S.
+///
+/// \copydetails ConvertToTH3I(const RHist<int> &hist)
+std::unique_ptr<TH3S> ConvertToTH3S(const RHist<short> &hist);
+
+/// Convert a two-dimensional histogram to TH3I.
+///
+/// If the RHistStats are tainted, for example after setting bin contents, the number of entries and the total sum of
+/// weights will be unset.
+///
+/// Throws an exception if the histogram has more than one dimension.
+///
+/// \param[in] hist the RHist to convert
+/// \return the converted TH3
+std::unique_ptr<TH3I> ConvertToTH3I(const RHist<int> &hist);
+
+/// Convert a two-dimensional histogram to TH3L.
+///
+/// \copydetails ConvertToTH3I(const RHist<int> &hist)
+std::unique_ptr<TH3L> ConvertToTH3L(const RHist<long> &hist);
+
+/// Convert a two-dimensional histogram to TH3L.
+///
+/// \copydetails ConvertToTH3I(const RHist<int> &hist)
+std::unique_ptr<TH3L> ConvertToTH3L(const RHist<long long> &hist);
+
+/// Convert a two-dimensional histogram to TH3F.
+///
+/// \copydetails ConvertToTH3I(const RHist<int> &hist)
+std::unique_ptr<TH3F> ConvertToTH3F(const RHist<float> &hist);
+
+/// Convert a two-dimensional histogram to TH3D.
+///
+/// \copydetails ConvertToTH3I(const RHist<int> &hist)
+std::unique_ptr<TH3D> ConvertToTH3D(const RHist<double> &hist);
+
+/// Convert a two-dimensional histogram to TH3D.
+///
+/// \copydetails ConvertToTH3I(const RHist<int> &hist)
+std::unique_ptr<TH3D> ConvertToTH3D(const RHist<RBinWithError> &hist);
+
+} // namespace Hist
+} // namespace Experimental
+} // namespace ROOT
+
+#endif

hist/histv7util/src/ConvertToTH3.cxx

@@ -0,0 +1,246 @@
+/// \file
+/// \warning This is part of the %ROOT 7 prototype! It will change without notice. It might trigger earthquakes.
+/// Feedback is welcome!
+
+#include <ROOT/Hist/ConversionUtils.hxx>
+#include <ROOT/Hist/ConvertToTH3.hxx>
+#include <ROOT/RBinIndex.hxx>
+#include <ROOT/RHist.hxx>
+#include <ROOT/RHistEngine.hxx>
+
+#include <TH3.h>
+
+#include <cassert>
+#include <cstddef>
+#include <memory>
+#include <stdexcept>
+#include <type_traits>
+#include <variant>
+#include <vector>
+
+using namespace ROOT::Experimental;
+
+namespace {
+template <typename Hist, typename T>
+std::unique_ptr<Hist> ConvertToTH3Impl(const RHistEngine<T> &engine)
+{
+   if (engine.GetNDimensions() != 3) {
+      throw std::invalid_argument("TH3 requires two dimensions");
+   }
+
+   auto ret = std::make_unique<Hist>();
+   ret->SetDirectory(nullptr);
+
+   const auto &axis0 = engine.GetAxes()[0];
+   ROOT::Experimental::Hist::Internal::ConvertAxis(*ret->GetXaxis(), axis0);
+   const auto &axis1 = engine.GetAxes()[1];
+   ROOT::Experimental::Hist::Internal::ConvertAxis(*ret->GetYaxis(), axis1);
+   const auto &axis2 = engine.GetAxes()[2];
+   ROOT::Experimental::Hist::Internal::ConvertAxis(*ret->GetZaxis(), axis2);
+   ret->SetBinsLength();
+
+   Double_t *sumw2 = nullptr;
+   auto copyBinContent = [&ret, &engine, &sumw2](Int_t i, RBinIndex index0, RBinIndex index1, RBinIndex index2) {
+      if constexpr (std::is_same_v<T, RBinWithError>) {
+         if (sumw2 == nullptr) {
+            ret->Sumw2();
+            sumw2 = ret->GetSumw2()->GetArray();
+         }
+         const RBinWithError &c = engine.GetBinContent(index0, index1, index2);
+         ret->GetArray()[i] = c.fSum;
+         sumw2[i] = c.fSum2;
+      } else {
+         (void)sumw2;
+         ret->GetArray()[i] = engine.GetBinContent(index0, index1, index2);
+      }
+   };
+
+   // Copy the bin contents, accounting for TH3 numbering conventions.
+   for (auto index0 : axis0.GetFullRange()) {
+      Int_t i0 = 0;
+      if (index0.IsUnderflow()) {
+         i0 = 0;
+      } else if (index0.IsOverflow()) {
+         i0 = axis0.GetNNormalBins() + 1;
+      } else {
+         assert(index0.IsNormal());
+         i0 = index0.GetIndex() + 1;
+      }
+      Int_t n0 = ret->GetXaxis()->GetNbins() + 2;
+
+      for (auto index1 : axis1.GetFullRange()) {
+         Int_t i1 = 0;
+         if (index1.IsUnderflow()) {
+            i1 = 0;
+         } else if (index1.IsOverflow()) {
+            i1 = axis1.GetNNormalBins() + 1;
+         } else {
+            assert(index1.IsNormal());
+            i1 = index1.GetIndex() + 1;
+         }
+         Int_t n1 = ret->GetYaxis()->GetNbins() + 2;
+
+         for (auto index2 : axis2.GetFullRange()) {
+            Int_t i2 = 0;
+            if (index2.IsUnderflow()) {
+               i2 = 0;
+            } else if (index2.IsOverflow()) {
+               i2 = axis2.GetNNormalBins() + 1;
+            } else {
+               assert(index2.IsNormal());
+               i2 = index2.GetIndex() + 1;
+            }
+
+            Int_t i = i0 + n0 * (i1 + n1 * i2);
+            copyBinContent(i, index0, index1, index2);
+         }
+      }
+   }
+
+   return ret;
+}
+
+template <typename Hist>
+void ConvertGlobalStatistics(Hist &h, const RHistStats &stats)
+{
+   if (stats.IsTainted()) {
+      return;
+   }
+
+   h.SetEntries(stats.GetNEntries());
+
+   Double_t hStats[11] = {
+      stats.GetSumW(),
+      stats.GetSumW2(),
+      0,
+      0,
+      0,
+      0,
+      // We do not have sumwxy
+      0,
+      0,
+      0,
+      // We do not have sumwxz
+      0,
+      // We do not have sumwyz
+      0,
+   };
+   if (stats.IsEnabled(0)) {
+      hStats[2] = stats.GetDimensionStats(0).fSumWX;
+      hStats[3] = stats.GetDimensionStats(0).fSumWX2;
+   }
+   if (stats.IsEnabled(1)) {
+      hStats[4] = stats.GetDimensionStats(1).fSumWX;
+      hStats[5] = stats.GetDimensionStats(1).fSumWX2;
+   }
+   if (stats.IsEnabled(2)) {
+      hStats[7] = stats.GetDimensionStats(2).fSumWX;
+      hStats[8] = stats.GetDimensionStats(2).fSumWX2;
+   }
+   h.PutStats(hStats);
+}
+} // namespace
+
+namespace ROOT {
+namespace Experimental {
+namespace Hist {
+
+std::unique_ptr<TH3C> ConvertToTH3C(const RHistEngine<char> &engine)
+{
+   return ConvertToTH3Impl<TH3C>(engine);
+}
+
+std::unique_ptr<TH3S> ConvertToTH3S(const RHistEngine<short> &engine)
+{
+   return ConvertToTH3Impl<TH3S>(engine);
+}
+
+std::unique_ptr<TH3I> ConvertToTH3I(const RHistEngine<int> &engine)
+{
+   return ConvertToTH3Impl<TH3I>(engine);
+}
+
+std::unique_ptr<TH3L> ConvertToTH3L(const RHistEngine<long> &engine)
+{
+   return ConvertToTH3Impl<TH3L>(engine);
+}
+
+std::unique_ptr<TH3L> ConvertToTH3L(const RHistEngine<long long> &engine)
+{
+   return ConvertToTH3Impl<TH3L>(engine);
+}
+
+std::unique_ptr<TH3F> ConvertToTH3F(const RHistEngine<float> &engine)
+{
+   return ConvertToTH3Impl<TH3F>(engine);
+}
+
+std::unique_ptr<TH3D> ConvertToTH3D(const RHistEngine<double> &engine)
+{
+   return ConvertToTH3Impl<TH3D>(engine);
+}
+
+std::unique_ptr<TH3D> ConvertToTH3D(const RHistEngine<RBinWithError> &engine)
+{
+   return ConvertToTH3Impl<TH3D>(engine);
+}
+
+std::unique_ptr<TH3C> ConvertToTH3C(const RHist<char> &hist)
+{
+   auto ret = ConvertToTH3C(hist.GetEngine());
+   ConvertGlobalStatistics(*ret, hist.GetStats());
+   return ret;
+}
+
+std::unique_ptr<TH3S> ConvertToTH3S(const RHist<short> &hist)
+{
+   auto ret = ConvertToTH3S(hist.GetEngine());
+   ConvertGlobalStatistics(*ret, hist.GetStats());
+   return ret;
+}
+
+std::unique_ptr<TH3I> ConvertToTH3I(const RHist<int> &hist)
+{
+   auto ret = ConvertToTH3I(hist.GetEngine());
+   ConvertGlobalStatistics(*ret, hist.GetStats());
+   return ret;
+}
+
+std::unique_ptr<TH3L> ConvertToTH3L(const RHist<long> &hist)
+{
+   auto ret = ConvertToTH3L(hist.GetEngine());
+   ConvertGlobalStatistics(*ret, hist.GetStats());
+   return ret;
+}
+
+std::unique_ptr<TH3L> ConvertToTH3L(const RHist<long long> &hist)
+{
+   auto ret = ConvertToTH3L(hist.GetEngine());
+   ConvertGlobalStatistics(*ret, hist.GetStats());
+   return ret;
+}
+
+std::unique_ptr<TH3F> ConvertToTH3F(const RHist<float> &hist)
+{
+   auto ret = ConvertToTH3F(hist.GetEngine());
+   ConvertGlobalStatistics(*ret, hist.GetStats());
+   return ret;
+}
+
+std::unique_ptr<TH3D> ConvertToTH3D(const RHist<double> &hist)
+{
+   auto ret = ConvertToTH3D(hist.GetEngine());
+   ConvertGlobalStatistics(*ret, hist.GetStats());
+   return ret;
+}
+
+std::unique_ptr<TH3D> ConvertToTH3D(const RHist<RBinWithError> &hist)
+{
+   auto ret = ConvertToTH3D(hist.GetEngine());
+   ConvertGlobalStatistics(*ret, hist.GetStats());
+   return ret;
+}
+
+} // namespace Hist
+} // namespace Experimental
+} // namespace ROOT

hist/histv7util/test/CMakeLists.txt

@@ -1,3 +1,4 @@
 ROOT_ADD_GTEST(hist_conversion_utils hist_conversion_utils.cxx LIBRARIES ROOTHist ROOTHistUtil)
 ROOT_ADD_GTEST(hist_convert_TH1 hist_convert_TH1.cxx LIBRARIES ROOTHist ROOTHistUtil)
 ROOT_ADD_GTEST(hist_convert_TH2 hist_convert_TH2.cxx LIBRARIES ROOTHist ROOTHistUtil)
+ROOT_ADD_GTEST(hist_convert_TH3 hist_convert_TH3.cxx LIBRARIES ROOTHist ROOTHistUtil)