Issue with RooGenericPdf::createIntegral in 2D

Hello,

I’ve built a RooGenericPdf and I want to find the normalized integral of it within certain bounds. I want to fill the bins of a TH2F with these values where the bin edges serve as the bounds of my pdf. I do this by specifying:

RooGeneric pdf("pdf","pdf", someformulaofXY, RooArgList(x,y));
RooAbsReal* iVar = pdf.createIntegral(RooArgSet(x,y),NormSet(RooArgSet(x,y)),Range("thisbin"));

// Enter Loop
for(int i=1; ...) {
for(int j=1; ...) {
   x.setRange("thisbin", binLoX, binHiX);
   y.setRange("thisbin", binLoY, binHiY);
   double NormalizedIntegralOverBin = iVar->getVal();
   hist2D->SetBinContent(i,j, NormalizedIntegralOverBin);
}
}

I use a loop to loop over every bin, set the appropriate “thisbin” range for x and y, and fill the bin with iVar->getVal(). The problem is that the value of “iVar->getVal()” doesnt change as I change the range [min, max] for x and y.

I’ve checked this result by modifying a copy of “rf110_normintegration.C” and noticed that everything works fine for a RooGaussian and a RooGenericPdf in 1D, but when I go to 2D only the RooGaussian seems to work. The integral of the 2D RooGenericPdf always returns the same value (the value I get from the the first time I set the range).

I’ve tried a number of things with no luck including the most obvious ones:

  • Create/Delete the integral variable in each loop after setting the ranges (didnt work)
  • Use “x.removeRange(“thisbin”)” at the end of each loop (also didnt work)
    The only solution I can come up with is to create a unique range in each iteration (something like “thisbin_i_j” where i,j are bin indices) and create the integral variable using this unique range. However this is very tedious and takes forever with the number of bins I’m trying to use. Not to mention the obvious memory leak storing all these ranges.

I’ve attached copies of the files which show this error (both 1D and 2D). Am I doing something wrong? Is there a better way to do this? Any help anyone can give me would be very much appreciated.

Thanks!
Josh

=== EDIT ===
I should note that I’m using Root v5.34.04, Roofit 3.55 on Ubuntu 12.10.
rf110_normintegration_2Dmod.C (2.61 KB)
rf110_normintegration_1Dmod.C (2.14 KB)

Update:

So I was curious about this issue, so I kept poking at it. I wanted to test if the fact that the RooGaussian was being integrated analytically while the RooGenericPdf was being integrated numerically had anything to do with the behavior I was getting. So I after defining the RooProdPdf of the two RooGaussian’s I forced it to intgrate the 2D Gaussian numerically using:

  RooGaussian gx("gx","gx",x,RooConst(-2),RooConst(3)) ;
  RooGaussian gy("gy","gy",y,RooConst(-2),RooConst(3)) ;
  RooProdPdf gxy("gxy","gxy", gx, gy);
  gxy.forceNumInt(kTRUE);

and now the 2D Gaussian has the same behavior as the RooGenericPdf. So it must have something to do with the pdf’s NUMERIC integration. I dont know if this helps at all, but I figured I would pass it along just in case.

In the mean time I’ve just been using a unique range name for each iteration of the loop. It consumes more memory though, so if there is a way to fix this issue I would love to hear it. Thanks again for any help you can give!

  • Josh

I’m going to continue documenting my progress on this issue as I go.

I noticed that what is actually being returned by RooGenericPdf::createIntegral is a RooRealIntegral “disguised” as a RooAbsReal. Also, by default RooRealIntegral caches all 2D evaluations (I assume for speed reasons). This had the negative impact that I was seeing. It is possible to set how many dimensions are required in your PDF before a RooRealIntegral starts caching the value via RooRealIntegral::setCacheAllNumeric. So, in order to fix my problem I changed my old code to include the following line at the beginning:

RooRealIntegral::setCacheAllNumeric(3);

Now, I’ve required PDF’s to have 3 dimensions before the integral starts caching its value. I’ve attached an updated version of the 2D script from above showing this improvement.

Great! So this solves my problems and no longer requires me to set a unique range in each iteration of the loop. I can set this up before the loop and just call it after I’ve set the appropriate ranges on x,y under the name “thisbin”. Unfortunately the PDF’s I’m actually working with are very complicated and using this method increases the computation time by over an order of magnitude (it would take days to get the final result). I’m assuming this is due to the normalization step because when I remove the “NormSet()” parameter in “createIntegral” It’s actually faster by almost an order of magnitude. Obviously this gives me an incorrect result so I cant do that and its back to square one.

I believe I have mashed together a method for now which is sufficient for the specific case I am currently working with, but it would be nice in the future to know of a general way of doing this for any case. So if anyone can come up with a solution to this problem I would be greatly appreciative.