Integrating a RooAbsPdf derived class

[aaronsw]

Hi, I would like to write a derived RooAbsPdf class, “A” inside an extended RooAddPdf, “B”. When “B” gets evaluated, it integrates class “A” over the variable range:

>>> B->getVal();
## A::evaluate() called
## A::evaluate() called
## A::evaluate() called
...

The issue is that this doesn’t seem to stop or converge. What I see is that after calling B->getVal(), the program gets stuck calling A->evaluate() over and over again. The shape of A is a bit uneven since it is taken directly from a histogram.

My question is:

1. Is there any way to make this work? Is there a reason that this (guessed) integral of A does not converge? Is it due to sharp edges in the histogram?
2. If #1 isn’t possible, is it possible to provide the integral of A myself? Such that whenever B wants an integral, it calls A::getAnalyticalIntegral() or a similar function? What functions need to be implemented?

Thanks for any help. I can’t easily make a simple stand-alone example of this problem, but if is necessary I can try.

[StephanH]

Hi @aaronsw,

@aaronsw

I can’t easily make a simple stand-alone example of this problem, but if is necessary I can try.

That doesn’t make things easier. It’s not clear what you are doing.

What kind of integral to you want? Indefinite / definite / CDF?
You can find examples in your ${ROOTSYS}/tutorials/roofit/rf110_normintegration.C or
https://root.cern.ch/root/html/tutorials/roofit/rf110_normintegration.C.html

• If it’s definite, it should be easy. No writing of new classes:

   RooRealVar x("x", "x", -10, 10);
   RooGaussian gx("gx", "gx", x, RooConst(-2), RooConst(3));

   // Create object representing integral over gx
   // which is used to calculate  gx_Norm[x] == gx / gx_Int[x]
   auto igx = gx.createIntegral(x);
   cout << "gx_Int[x] = " << igx->getVal() << endl;

If you use RooFit::Range(<name of the range>), you can also define the integration range.

• CDFs are also easy. From the same example:

   // Create the cumulative distribution function of gx
   // i.e. calculate Int[-10,x] gx(x') dx'
   auto gx_cdf = gx.createCdf(x);

[aaronsw]

Hi @StephanH,

Thanks very much for your replies - I’ve spent a long time on these problems and your help has been invaluable. I put together an example of the problem: demo.zip (131.6 KB). The summary is that I want to make a PDF who’s shape is drawn from a TH1. I know there are other ways to do this in RooFit, but for the final model these won’t work.

Here there are three files:

• data.root: contains a histogram TH1
• simpleHist.h: a RooAbsPdf derived class called simpleHist that stores the TH1 and returns the value of the TH1 as the shape of the PDF (in evaluate()):

        int binIndex = hist->FindBin(x);
        double ret   = hist->GetBinContent(binIndex)/hist->GetBinWidth(binIndex);
        return ret;

• demo.py which makes two PDF’s: testPdf and exp. It then adds these PDF’s together in an extended RooAbsPdf called addedPdfs. Finally it calls addedPdf.getVal(), which triggers the never ending loop of evaluations of testPdf

To see the problem, run demo.py. When it reaches print addedPdfs.getVal() the infinite loop starts. Since the evaluate() of testPdf keeps getting called, I guess this is a loop trying to integrate testPdf (but I am unsure). (Interestingly, if the range of the extended PDF is the same as the range of the variable, this infinite loop doesn’t take place. Maybe this is because RooAddPdf doesn’t bother with the integral…) If anything is unclear or doesn’t work in this example, please let me know. This is the smallest implementation I could make that has this problem.

If you can think of any way to avoid this infinite loop, I think that would solve this problem. Thanks very much and apologies for the messy nature of this problem.

[StephanH]

Hi @aaronsw,

I just looked at it for 2 minutes, but here are some observations:

• The PDF is sometimes negative. That’s impossible for a PDF, and hence, RooFit forces the value to zero when you return a negative value. That means that Minuit will see a vanishing gradient, and won’t know in which direction to go next, because all the values that are returned are forced to 0.

[#0] ERROR:Eval -- RooAbsReal::logEvalError(f) evaluation error, 
 origin       : simpleHist::f[ x=x ]
 message      : p.d.f value is less than zero (-17.797049), forcing value to zero

• The PDF does not have a function for an analytical integral. Therefore, RooFit will start scanning the PDF to compute it’s integral numerically. This should then be used to normalise the PDF.

These two could already be the cause of the problems. If this doesn’t help, I could try to invest a bit more time to figure it out with you.

[aaronsw]

Hi @StephanH,

Thanks for the help. I’m not doing a fit here, so I think the negative values are okay. After a lot of experimentation, I reached the conclusion that a possible cause is the fact that the PDF based on a histogram is discontinuous. I made another stand-alone test of this: demo.zip (132.2 KB). By interpolating between the bins I was able to avoid the problem, but not really solve it.

Thanks

[StephanH]

Hi,

it makes sense if the steps (e.g. for adaptive numeric integration) are too small to make a difference. It will just see constant results everywhere.