From 4b05e75eea75993890cccd0f94dcef1491be0599 Mon Sep 17 00:00:00 2001
From: Jonas Rembser <jonas.rembser@cern.ch>
Date: Fri, 15 Jul 2022 13:55:54 +0200
Subject: [PATCH] [RF] Fixup to parameter index calculation in ParamHistFunc

Last year, with commit 3657e7c, the parameter index calculation was
changed to be on the fly instead of using a look-up map, which is much
faster.

However, the implemented formula was not correct for two or three
dimensions, which is fixed by this commit.

To make sure that the index computation is correct this time, the new
code was tested in this code snippet with various inputs:

```C++
void runTest(int nx = 42, int ny = 42, int nz = 42) {
  const int nxy = nx * ny;
  const int nyz = ny * nz;

  for (int i = 0; i < nx; ++i) {
    for (int j = 0; j < ny; ++j) {
      for (int k = 0; k < nz; ++k) {
        const int index = k + j * nz + i * ny * nz;
        const int gammaIndex = i + j * nx + k * nx * ny;

        const int i2 = index / nyz;
        const int tmp = index % nyz;
        const int j2 = tmp / nz;
        const int k2 = tmp % nz;

        const int gammaIndex2 = i2 + j2 * nx + k2 * nxy;

        if (gammaIndex2 != gammaIndex) {
          std::cout << "The unraveled indices were not correct!"
                    << std::endl;
          return;
        }
      }
    }
  }
}
```

Needs to be backported to the 6.26 branch to get into the 6.26.06 patch
release.

This commit the following problem reported on the forum:
https://root-forum.cern.ch/t/cpycppyy-segfault-on-mac-m1/50822
---
 roofit/histfactory/src/ParamHistFunc.cxx | 10 +++++++---
 1 file changed, 7 insertions(+), 3 deletions(-)

diff --git a/roofit/histfactory/src/ParamHistFunc.cxx b/roofit/histfactory/src/ParamHistFunc.cxx
index c3eb8bc7090c3..2aa87eef7245f 100644
--- a/roofit/histfactory/src/ParamHistFunc.cxx
+++ b/roofit/histfactory/src/ParamHistFunc.cxx
@@ -210,9 +210,13 @@ RooAbsReal& ParamHistFunc::getParameter( Int_t index ) const {
     _numBinsPerDim = getNumBinsPerDim(_dataVars);
   }
 
-  int i = index / n.yz ;
-  int j = (index % n.y) / n.z;
-  int k = index % (n.yz);
+  // Unravel the index to 3D coordinates. We can't use the index directly,
+  // because in the parameter set the dimensions are ordered in reverse order
+  // compared to the RooDataHist (for historical reasons).
+  const int i = index / n.yz;
+  const int tmp = index % n.yz;
+  const int j = tmp / n.z;
+  const int k = tmp % n.z;
 
   return static_cast<RooAbsReal&>(_paramSet[i + j * n.x + k * n.xy]);
 }