RooAbsPDF analytical integration doesn't seem to be used

I was having some issues normalizing some custom PDFs in roofit, and since the PDF should actually have an integral that never changes, I figured I would try to implement the integral myself. Unfortunately, my implementation never seems to be used. Here is my attempt at a dummy implementation, to see if I could get it to work:

from the pdf header file “RooSignalPDF.h”

               Int_t getAnalyticalIntegral(const RooArgSet &integSet,
                                           RooArgSet &anaIntSet,
                                           const char* rangeName=0) const;
               Double_t analyticalIntegral( Int_t code,
                                            const char* rangeName=0) const;

from the pdf source file “RooSignalPDF.cxx”


Int_t RooSignalPDF::getAnalyticalIntegral(const RooArgSet& integSet,
                                          RooArgSet& anaIntSet,
                                          const char* rangeName) const{
  cout << "**************do i even get here?************" << endl;
  return 1;  // always use this
}

// CAUTION: this is only valid as long as we are simply reading our model
// from pre-defined files that cannot change mid-run!
Double_t RooSignalPDF::analyticalIntegral(Int_t code,
                                          const char* rangeName) const{

  cout << "using analytical Integral" << endl;
  return 2;
}

and my attempt to call it:

  cout << "ci: " << ws->var("ci")->getVal() << endl;
  cout << "evts_per_ci: " << ws->var("evts_per_ci")->getVal() << endl;
  cout << "S1: " << ws->var("S1")->getVal() << endl;
  cout << "S2: " << ws->var("S2")->getVal() << endl;
  cout << "r: " << ws->var("r")->getVal() << endl;
  cout << "phi: " << ws->var("phi")->getVal() << endl;
  cout << "drift: " << ws->var("drift")->getVal() << endl;
  cout << endl;

  TString sigModName =
    TString::Format("sigPopTB%dZS%d", 1, 1);
  cout << "signalModel: " << ws->pdf(sigModName.Data())->getVal() << endl;
  cout << "signalModel norm: " << ws->pdf(sigModName.Data())->getNorm() << endl;

which gives the output

ci: 3.115e-09
evts_per_ci: 2.40713e+08
S1: 10
S2: 747.5
r: 9
phi: 0.00124981
drift: 57.0255

signalModel: 1.6928e-07
signalModel norm: 1

I would have expected the cout statements to show up somewhere in the above output, but they don’t.

My understanding from https://root.cern.ch/doc/v608/classRooAbsPdf.html is that whenever a normalization is done (or indeed even a “getVal()”) getAnalyticalIntegral() is called, and unless it returns 0, it calls analyticalIntegral() and takes its output as the result for the normalization integral.

Is this untrue, or did I do something wrong in my implementation?

Thanks so much!

Shaun Alsum

Hi @SAlsum,

Our expert @moneta will have a look when he is back on Monday. Apologies for the wait,

Enric

I did a little more digging and found the following in the getVal() documentation:
“Return current value, normalizated by integrating over the observables in ‘nset’. If ‘nset’ is 0, the unnormalized value. is returned.”
(which the . after value looks like a typo in the documention, just as an aside).

Which I imagine is why the normalization is not being called in this example. However, this leads me to another question, which is how to call it explicitly with an arg set. Simply trying to do it the naive way

ws->pdf(sigModName.Data())->getVal(
        RooArgSet(ws->var("S1"),
                             ws->var("S2"),
                             ws->var("r"),
                             ws->var("phi"),
                             ws->var("drift"))

leads to an error “/usr/common/usg/software/ROOT/6.08.00/include/root/RooArgSet.h:53:3: note: no known conversion for argument 1 from âRooReal
Var*â to âconst RooAbsArg&â”

and when I try to cast it myself via (RooAbsArg*)(ws->var(“S1”)), i get an invalid cast error.

Thank you so much for your help.

your code is not correct. RooWorkspace::var("…") returns a RooRealVar pointer. Here is the correct code:

ws->pdf(sigModName.Data())->getVal(
        RooArgSet(*ws->var("S1"),
                             *ws->var("S2"),
                             *ws->var("r"),
                             *ws->var("phi"),
                             *ws->var("drift"))

Lorenzo

This topic was automatically closed 14 days after the last reply. New replies are no longer allowed.