[hist] Implement axis slicing

hist/histv7/inc/ROOT/RAxisVariant.hxx

@@ -127,6 +127,26 @@ public:
       }
    }
 
+   /// Slice this axis according to the specification.
+   ///
+   /// Axes throw exceptions if the slicing cannot be performed. For example, the rebin operation must divide the
+   /// number of normal bins for RRegularAxis and RVariableBinAxis, while RCategoricalAxis cannot be sliced at all.
+   ///
+   /// \param[in] sliceSpec the slice specification
+   /// \return the sliced axis
+   RAxisVariant Slice(const RSliceSpec &sliceSpec) const
+   {
+      if (auto *regular = GetRegularAxis()) {
+         return regular->Slice(sliceSpec);
+      } else if (auto *variable = GetVariableBinAxis()) {
+         return variable->Slice(sliceSpec);
+      } else if (auto *categorical = GetCategoricalAxis()) {
+         return categorical->Slice(sliceSpec);
+      } else {
+         throw std::logic_error("unimplemented axis type"); // GCOVR_EXCL_LINE
+      }
+   }
+
    friend bool operator==(const RAxisVariant &lhs, const RAxisVariant &rhs) { return lhs.fVariant == rhs.fVariant; }
 
    /// %ROOT Streamer function to throw when trying to store an object of this class.

hist/histv7/inc/ROOT/RCategoricalAxis.hxx

@@ -8,6 +8,7 @@
 #include "RBinIndex.hxx"
 #include "RBinIndexRange.hxx"
 #include "RLinearizedIndex.hxx"
+#include "RSliceSpec.hxx"
 
 #include <cassert>
 #include <cstddef>
@@ -167,6 +168,31 @@ public:
                                 : GetNormalRange();
    }
 
+   /// Slice this axis according to the specification.
+   ///
+   /// A categorical axis cannot be sliced. The method will throw if a specification other than the default slice
+   /// operation is passed.
+   ///
+   /// \param[in] sliceSpec the slice specification
+   /// \return the sliced / copied axis
+   RCategoricalAxis Slice(const RSliceSpec &sliceSpec) const
+   {
+      if (sliceSpec.GetOperationSum() != nullptr) {
+         throw std::runtime_error("sum operation makes dimension disappear");
+      }
+
+      if (!sliceSpec.GetRange().IsInvalid()) {
+         throw std::runtime_error("slicing of RCategoricalAxis not implemented");
+      }
+      if (sliceSpec.GetOperationRebin() != nullptr) {
+         throw std::runtime_error("cannot rebin RCategoricalAxis");
+      }
+
+      // The sliced axis always has flow bins enabled, for symmetry with other axis types.
+      bool enableOverflowBin = true;
+      return RCategoricalAxis(fCategories, enableOverflowBin);
+   }
+
    /// %ROOT Streamer function to throw when trying to store an object of this class.
    void Streamer(TBuffer &) { throw std::runtime_error("unable to store RCategoricalAxis"); }
 };

hist/histv7/inc/ROOT/RRegularAxis.hxx

@@ -8,6 +8,7 @@
 #include "RBinIndex.hxx"
 #include "RBinIndexRange.hxx"
 #include "RLinearizedIndex.hxx"
+#include "RSliceSpec.hxx"
 
 #include <cassert>
 #include <cmath>
@@ -83,6 +84,24 @@ public:
    double GetHigh() const { return fHigh; }
    bool HasFlowBins() const { return fEnableFlowBins; }
 
+   /// Compute the low edge of a bin.
+   ///
+   /// \param[in] bin the index, must be \f$< fNNormalBins\f$
+   double ComputeLowEdge(std::uint64_t bin) const
+   {
+      assert(bin < fNNormalBins);
+      return fLow + (fHigh - fLow) * bin / fNNormalBins;
+   }
+
+   /// Compute the high edge of a bin.
+   ///
+   /// \param[in] bin the index, must be \f$< fNNormalBins\f$
+   double ComputeHighEdge(std::uint64_t bin) const
+   {
+      assert(bin < fNNormalBins);
+      return fLow + (fHigh - fLow) * (bin + 1) / fNNormalBins;
+   }
+
    friend bool operator==(const RRegularAxis &lhs, const RRegularAxis &rhs)
    {
       return lhs.fNNormalBins == rhs.fNNormalBins && lhs.fLow == rhs.fLow && lhs.fHigh == rhs.fHigh &&
@@ -194,6 +213,59 @@ public:
                              : GetNormalRange();
    }
 
+   /// Slice this axis according to the specification.
+   ///
+   /// Throws an exception if the axis cannot be sliced:
+   ///  * A sum operation makes the dimension disappear.
+   ///  * The rebin operation must divide the number of normal bins.
+   ///
+   /// \param[in] sliceSpec the slice specification
+   /// \return the sliced axis, with enabled underflow and overflow bins
+   RRegularAxis Slice(const RSliceSpec &sliceSpec) const
+   {
+      if (sliceSpec.GetOperationSum() != nullptr) {
+         throw std::runtime_error("sum operation makes dimension disappear");
+      }
+
+      // Figure out the properties of the sliced axis.
+      std::uint64_t nNormalBins = fNNormalBins;
+      double low = fLow;
+      double high = fHigh;
+
+      const auto &range = sliceSpec.GetRange();
+      if (!range.IsInvalid()) {
+         std::uint64_t begin = 0;
+         std::uint64_t end = nNormalBins;
+         if (range.GetBegin().IsNormal()) {
+            begin = range.GetBegin().GetIndex();
+            // Only compute a new lower end of the axis interval if needed.
+            if (begin > 0) {
+               low = ComputeLowEdge(begin);
+            }
+         }
+         if (range.GetEnd().IsNormal()) {
+            end = range.GetEnd().GetIndex();
+            // Only compute a new upper end of the axis interval if needed, to avoid floating-point inaccuracies.
+            if (end < nNormalBins) {
+               high = ComputeHighEdge(end - 1);
+            }
+         }
+         nNormalBins = end - begin;
+      }
+
+      if (auto *opRebin = sliceSpec.GetOperationRebin()) {
+         if (nNormalBins % opRebin->GetNGroup() != 0) {
+            throw std::runtime_error("rebin operation does not divide number of normal bins");
+         }
+         nNormalBins /= opRebin->GetNGroup();
+      }
+
+      // The sliced axis always has flow bins enabled to preserve all entries. This is the least confusing for users,
+      // even if not always strictly necessary.
+      bool enableFlowBins = true;
+      return RRegularAxis(nNormalBins, {low, high}, enableFlowBins);
+   }
+
    /// %ROOT Streamer function to throw when trying to store an object of this class.
    void Streamer(TBuffer &) { throw std::runtime_error("unable to store RRegularAxis"); }
 };

hist/histv7/inc/ROOT/RVariableBinAxis.hxx

@@ -8,6 +8,7 @@
 #include "RBinIndex.hxx"
 #include "RBinIndexRange.hxx"
 #include "RLinearizedIndex.hxx"
+#include "RSliceSpec.hxx"
 
 #include <cassert>
 #include <cmath>
@@ -196,6 +197,57 @@ public:
                              : GetNormalRange();
    }
 
+   /// Slice this axis according to the specification.
+   ///
+   /// Throws an exception if the axis cannot be sliced:
+   ///  * A sum operation makes the dimension disappear.
+   ///  * The rebin operation must divide the number of normal bins.
+   ///
+   /// \param[in] sliceSpec the slice specification
+   /// \return the sliced axis, with enabled underflow and overflow bins
+   RVariableBinAxis Slice(const RSliceSpec &sliceSpec) const
+   {
+      if (sliceSpec.GetOperationSum() != nullptr) {
+         throw std::runtime_error("sum operation makes dimension disappear");
+      }
+
+      // Figure out the properties of the sliced axis.
+      std::size_t nNormalBins = fBinEdges.size() - 1;
+      std::size_t begin = 0;
+      std::size_t end = nNormalBins;
+
+      const auto &range = sliceSpec.GetRange();
+      if (!range.IsInvalid()) {
+         if (range.GetBegin().IsNormal()) {
+            begin = range.GetBegin().GetIndex();
+         }
+         if (range.GetEnd().IsNormal()) {
+            end = range.GetEnd().GetIndex();
+         }
+         nNormalBins = end - begin;
+      }
+
+      std::uint64_t nGroup = 1;
+      if (auto *opRebin = sliceSpec.GetOperationRebin()) {
+         nGroup = opRebin->GetNGroup();
+         if (nNormalBins % nGroup != 0) {
+            throw std::runtime_error("rebin operation does not divide number of normal bins");
+         }
+      }
+
+      // Prepare the bin edges.
+      std::vector<double> binEdges;
+      for (std::size_t bin = begin; bin < end; bin += nGroup) {
+         binEdges.push_back(fBinEdges[bin]);
+      }
+      binEdges.push_back(fBinEdges[end]);
+
+      // The sliced axis always has flow bins enabled to preserve all entries. This is the least confusing for users,
+      // even if not always strictly necessary.
+      bool enableFlowBins = true;
+      return RVariableBinAxis(std::move(binEdges), enableFlowBins);
+   }
+
    /// %ROOT Streamer function to throw when trying to store an object of this class.
    void Streamer(TBuffer &) { throw std::runtime_error("unable to store RVariableBinAxis"); }
 };

hist/histv7/test/hist_axis_variant.cxx

@@ -23,6 +23,11 @@ TEST(RAxisVariant, RegularAxis)
       EXPECT_EQ(std::distance(normal12.begin(), normal12.end()), 1);
       const auto full = axis.GetFullRange();
       EXPECT_EQ(std::distance(full.begin(), full.end()), Bins + 2);
+
+      const auto slicedAxis = axis.Slice(RSliceSpec{});
+      EXPECT_EQ(slicedAxis.GetVariant().index(), 0);
+      EXPECT_TRUE(slicedAxis.GetRegularAxis() != nullptr);
+      EXPECT_EQ(slicedAxis.GetNNormalBins(), Bins);
    }
 
    {
@@ -36,6 +41,11 @@ TEST(RAxisVariant, RegularAxis)
       EXPECT_EQ(std::distance(normal12.begin(), normal12.end()), 1);
       const auto full = axis.GetFullRange();
       EXPECT_EQ(std::distance(full.begin(), full.end()), Bins);
+
+      const auto slicedAxis = axis.Slice(RSliceSpec{});
+      EXPECT_EQ(slicedAxis.GetVariant().index(), 0);
+      EXPECT_TRUE(slicedAxis.GetRegularAxis() != nullptr);
+      EXPECT_EQ(slicedAxis.GetNNormalBins(), Bins);
    }
 }
 
@@ -64,6 +74,11 @@ TEST(RAxisVariant, VariableBinAxis)
       EXPECT_EQ(std::distance(normal12.begin(), normal12.end()), 1);
       const auto full = axis.GetFullRange();
       EXPECT_EQ(std::distance(full.begin(), full.end()), Bins + 2);
+
+      const auto slicedAxis = axis.Slice(RSliceSpec{});
+      EXPECT_EQ(slicedAxis.GetVariant().index(), 1);
+      EXPECT_TRUE(slicedAxis.GetVariableBinAxis() != nullptr);
+      EXPECT_EQ(slicedAxis.GetNNormalBins(), Bins);
    }
 
    {
@@ -77,6 +92,11 @@ TEST(RAxisVariant, VariableBinAxis)
       EXPECT_EQ(std::distance(normal12.begin(), normal12.end()), 1);
       const auto full = axis.GetFullRange();
       EXPECT_EQ(std::distance(full.begin(), full.end()), Bins);
+
+      const auto slicedAxis = axis.Slice(RSliceSpec{});
+      EXPECT_EQ(slicedAxis.GetVariant().index(), 1);
+      EXPECT_TRUE(slicedAxis.GetVariableBinAxis() != nullptr);
+      EXPECT_EQ(slicedAxis.GetNNormalBins(), Bins);
    }
 }
 
@@ -100,6 +120,11 @@ TEST(RAxisVariant, CategoricalAxis)
       EXPECT_EQ(std::distance(normal12.begin(), normal12.end()), 1);
       const auto full = axis.GetFullRange();
       EXPECT_EQ(std::distance(full.begin(), full.end()), 4);
+
+      const auto slicedAxis = axis.Slice(RSliceSpec{});
+      EXPECT_EQ(slicedAxis.GetVariant().index(), 2);
+      EXPECT_TRUE(slicedAxis.GetCategoricalAxis() != nullptr);
+      EXPECT_EQ(slicedAxis.GetNNormalBins(), 3);
    }
 
    {
@@ -113,5 +138,10 @@ TEST(RAxisVariant, CategoricalAxis)
       EXPECT_EQ(std::distance(normal12.begin(), normal12.end()), 1);
       const auto full = axis.GetFullRange();
       EXPECT_EQ(std::distance(full.begin(), full.end()), 3);
+
+      const auto slicedAxis = axis.Slice(RSliceSpec{});
+      EXPECT_EQ(slicedAxis.GetVariant().index(), 2);
+      EXPECT_TRUE(slicedAxis.GetCategoricalAxis() != nullptr);
+      EXPECT_EQ(slicedAxis.GetNNormalBins(), 3);
    }
 }

hist/histv7/test/hist_categorical.cxx

@@ -198,3 +198,35 @@ TEST(RCategoricalAxis, GetFullRange)
       EXPECT_EQ(std::distance(full.begin(), full.end()), categories.size());
    }
 }
+
+static void Test_RCategoricalAxis_Slice(bool enableOverflowBin)
+{
+   const std::vector<std::string> categories = {"a", "b", "c"};
+   const RCategoricalAxis origAxis(categories, enableOverflowBin);
+   ASSERT_EQ(origAxis.HasOverflowBin(), enableOverflowBin);
+
+   {
+      // The only allowed "slicing" is to keep the entire axis.
+      const RSliceSpec full;
+      const auto axis = origAxis.Slice(full);
+      EXPECT_EQ(axis.GetNNormalBins(), 3);
+      EXPECT_EQ(axis.GetCategories(), categories);
+      EXPECT_TRUE(axis.HasOverflowBin());
+   }
+
+   // Slicing and other operations are not allowed / implemented.
+   EXPECT_THROW(origAxis.Slice(origAxis.GetFullRange()), std::runtime_error);
+   EXPECT_THROW(origAxis.Slice(origAxis.GetNormalRange()), std::runtime_error);
+   EXPECT_THROW(origAxis.Slice(RSliceSpec::ROperationRebin(2)), std::runtime_error);
+   EXPECT_THROW(origAxis.Slice(RSliceSpec::ROperationSum{}), std::runtime_error);
+}
+
+TEST(RCategoricalAxis, Slice)
+{
+   Test_RCategoricalAxis_Slice(true);
+}
+
+TEST(RCategoricalAxis, SliceNoOverflowBin)
+{
+   Test_RCategoricalAxis_Slice(false);
+}