Interpolation between two 2D PDFs using RooMomentMorphFuncND fails because "integrating a RooAbsRealLValue is not allowed"

ROOT Version: 6.32.02
Built for linuxx8664gcc on Sep 18 2024, 20:01:03
From heads/master@tags/v6-32-02

Hello ROOT experts. I have a bunch of TH2 histograms representing signal mass distributions binned in the space of two invariant masses (let’s say m_X and m_Y). For a given m_Y, I would like to generate a new signal shape at X mass m_x by interpolating between the shapes of two other signals at X masses m_{x1} and m_{x2}, for m_{x1} < m_x < m_{x2}. I found several similar threads on the ROOT forum which make use of moment morphing:

In particular, I verified that I could reproduce the results in the example code from the first link above.

I’ve written a similar function based on the linked code and attached it in the file Interpolator.py. This file opens two (attached) ROOT files and gets the TH2 histograms for signals of m_Y=1800 and m_{x1} = 125 and m_{x2} = 175. The code then creates a RooDataHist for each distribution which is used to make a RooHistPdf. I then create a RooMomentMorphFuncND.Grid2 object to which I add the two PDFs and then attempt to perform linear interpolation between them at the given mass point (150, 1800) using RooMomentMorphFuncND.Linear. However, the code fails with the error:

# DEBUG----------------------------------------------------------------------------------
RooRealVar::mPhi = 310  L(60 - 560) 
 <class cppyy.gbl.RooRealVar at 0x5a1c3ee41dd0>
TH1:        68.97920879694776
RDH:        68.97920879694763
RDH (true): 68979.20879694764
RHP:        0.0
TH1:        64.96953476283132
RDH:        64.96953476283137
RDH (true): 64969.534762831376
RHP:        5.955384569195107e-06
# DEBUG----------------------------------------------------------------------------------
[#0] ERROR:InputArguments -- Attempting to integrate the RooRealVar "mPhi", but integrating a RooAbsRealLValue is not allowed!
Traceback (most recent call last):
  File "/home/amitav/work/CMS/Tprime/Interpolation/Interpolator.py", line 81, in <module>
    Interpolate()
  File "/home/amitav/work/CMS/Tprime/Interpolation/Interpolator.py", line 69, in Interpolate
    pdf.plotOn(framex, RooFit.LineColor(ROOT.kBlue))
  File "/home/amitav/miniconda3/envs/ROOT/lib/python3.12/site-packages/ROOT/_pythonization/_roofit/_rooabspdf.py", line 79, in plotOn
    return self._plotOn(*args, **kwargs)
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
cppyy.gbl.std.runtime_error: RooPlot* RooAbsPdf::plotOn(RooPlot* frame, const RooCmdArg& arg1 = {}, const RooCmdArg& arg2 = {}, const RooCmdArg& arg3 = {}, const RooCmdArg& arg4 = {}, const RooCmdArg& arg5 = {}, const RooCmdArg& arg6 = {}, const RooCmdArg& arg7 = {}, const RooCmdArg& arg8 = {}, const RooCmdArg& arg9 = {}, const RooCmdArg& arg10 = {}) =>
    runtime_error: Attempting to integrate the RooRealVar "mPhi", but integrating a RooAbsRealLValue is not allowed!

I understand what the error message is saying, but it’s not clear to me why I see this failure mode when the same RooRealVar objects worked in the example code. Am I making a mistake in my code somewhere, perhaps where I convert the RooDataHist to a RooHistPdf? The RooHistPdf.analyticalIntegral(0) call returning 0 doesn’t seem correct to me.

Any help would be much appreciated! I’ve attached a minimal reproducible example to this post as well. Thanks in advance :slight_smile:

Interpolator.py (2.5 KB)
THselection_TprimeB-1800-125_18.root (630.4 KB)
THselection_TprimeB-1800-175_18.root (746.7 KB)

The code is attached in the file Interpolator.py, but included here as well:

import ROOT
from ROOT import TFile, RooRealVar, RooArgList, RooDataHist, RooHistPdf, RooFit, RooArgSet, RooMomentMorphFuncND, RooBinning, RooWrapperPdf

# X = mPhi
xMin = 60
xMax = 560
nX   = 50

# Y = mT
yMin = 800
yMax = 3500
nY   = 27

def create_RDH_PDF(name, h):
    xvar = RooRealVar('mPhi','mPhi',xMin,xMax)
    yvar = RooRealVar('mT','mT',yMin,yMax)
    RAL_myvars = RooArgList(xvar,yvar)
    RAS_myvars = RooArgSet(RAL_myvars)
    RDH = RooDataHist('RDH','RDH',RAL_myvars,h)
    RHP = RooHistPdf('RHP','RHP',RAS_myvars,RDH)

    # DEBUG
    print('TH1:        ' + str(h.Integral()))
    print('RDH:        ' + str(RDH.sum(False)))
    print('RDH (true): ' + str(RDH.sum(True)))
    print('RHP:        ' + str(RHP.analyticalIntegral(0)))

    return RDH, RHP

def getBin(binning, val):
    if (val > binning.highBound()):
        return binning.numBins()
    else:
        return binning.binNumber(val)

def Interpolate():
    # Observables
    x = RooRealVar('mPhi','mPhi',xMin,xMax)
    y = RooRealVar('mT','mT',yMin,yMax)
    print(x, type(x))
    # Interpolated observables
    ix = RooRealVar('i_mPhi','i_mPhi',xMin,xMax)
    iy = RooRealVar('i_mT','i_mT',yMin,yMax)
    # Interpolated binning 
    RBx = RooBinning(nX,xMin,xMax)
    RBy = RooBinning(nY,yMin,yMax)
    grid = RooMomentMorphFuncND.Grid2(RBx,RBy)
    # Get the PDFs for the two masses we want to interpolate b/w
    f1 = TFile.Open('THselection_TprimeB-1800-125_18.root')
    f2 = TFile.Open('THselection_TprimeB-1800-175_18.root')
    h1 = f1.Get('MHvsMTH_SR_pass__nominal')
    h2 = f2.Get('MHvsMTH_SR_pass__nominal')
    RDH1, RHP1 = create_RDH_PDF('1',h1)
    RDH2, RHP2 = create_RDH_PDF('2',h2)
    # Add the PDFs to the grid at their X,Y location
    grid.addPdf(RHP1, getBin(RBx, 1800), getBin(RBy, 125))
    grid.addPdf(RHP2, getBin(RBx, 1800), getBin(RBy, 175))
    # Morph 
    morph = RooMomentMorphFuncND('morph', 'morph', RooArgList(ix,iy),RooArgList(x,y), grid, RooMomentMorphFuncND.Linear)
    morph.setPdfMode()
    pdf = RooWrapperPdf('morph_pdf', 'morph_pdf', morph, True)
    # Get the frames for plotting 
    framex = x.frame()
    framey = y.frame()
    # Generate at a given point 
    ix.setVal(1800)
    ix.setVal(150)
    # Plot
    pdf.plotOn(framex, RooFit.LineColor(ROOT.kBlue))
    pdf.plotOn(framey, RooFit.LineColor(ROOT.kBlue))
    c = ROOT.TCanvas('c','c')
    c.Divide(2,1)
    c.cd(1)
    framex.Draw()
    c.cd(2)
    framey.Draw()
    framey.Draw('same')
    c.Print('out.pdf')


Interpolate()

Hello,

this is a good question for @jonas (who is very busy this week).
In the mean time, you could try to use pdf.Print("T") to print the tree of PDFs/functions being evaluated. Maybe that helps to understand why it thinks it needs to integrate a variable.

RooHistPdf.analyticalIntegral(0) might just mean that there is no analytical integral for this PDF, and that it needs to be integrated numerically (and RooFit seems to be trying that).

Hi Stephan,

Thanks very much for the quick reply. I think the error was somewhere in the creation of the RooDataHists, but I’m not sure where exactly. I worked around it by fitting the two existing signals (125, 1800) and (175, 1800) with a product of RooCrystalBall functions along the two axes and saving the resulting RooProdPdf to a RooWorkspace in a ROOT file. This works well, as evidenced by the fit result plot of the PDFs:

However, when I attempt to perform the interpolation for, e.g, a signal at (150, 1800), it seems only to copy the signal at (125, 1800) instead of interpolating between mPhi=125 and mPhi=175.

There are some messages from RooFit printed to the terminal (the error messages I believe can be ignored, as per here):

------------------------------------------------------------------------------------------
Interpolating signal of mPhi = 150, between existing signals at mPhi = [125,175]
------------------------------------------------------------------------------------------
[#1] INFO:NumericIntegration -- RooRealIntegral::init(model_MOMENT_1_mPhi_product_Int[mPhi,mT]_Norm[mPhi,mT]) using numeric integrator RooAdaptiveIntegratorND to calculate Int(mPhi,mT)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(model_MOMENT_1_mT_product_Int[mPhi,mT]_Norm[mPhi,mT]) using numeric integrator RooAdaptiveIntegratorND to calculate Int(mT,mPhi)
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name model_MOMENT_2C_mPhi is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name model_MOMENT_2C_mT is already in this set
[#1] INFO:NumericIntegration -- RooRealIntegral::init(model_MOMENT_2C_mPhi_product_Int[mPhi,mT]_Norm[mPhi,mT]) using numeric integrator RooAdaptiveIntegratorND to calculate Int(mPhi,mT)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(model_MOMENT_2C_mT_product_Int[mPhi,mT]_Norm[mPhi,mT]) using numeric integrator RooAdaptiveIntegratorND to calculate Int(mT,mPhi)
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name model_MOMENT_2C_mPhi is already in this set
[#0] ERROR:InputArguments -- RooArgSet::checkForDup: ERROR argument with name model_MOMENT_2C_mT is already in this set
Info in <TCanvas::Print>: pdf file TEST.pdf has been created

and, interestingly, when I run Print("T") on the morphed PDF, I get the error

[#0] ERROR:InputArguments -- RooLinearVar::RooLinearVar(morph_func_transVar_0_0): ERROR, slope(morph_func_slope_0_0) and offset(morph_func_offset_0_0) may not depend on variable(mPhi)

which crashes the script. I’m not too sure how to interpret this, but since I’m trying to perform the interpolation along the variable mPhi, the fact that the morphed PDF may not depend on this variable is troubling!

I know this is slightly off-topic from my main post, and that the RooMomentMorphFuncND.Grid2 and RooMomentMorphFuncND classes in general are new and not very well-documented, but if you or perhaps Jonas (absolutely no rush) have any insights, I would greatly appreciate it.

I’ve attached the files needed to run the reproducible example, and the script can be run via

python InterpolateShapes.py -m 150 -m1 125 -m2 175 -mT 1800 -y 17

1800-125_17_ws.root (9.7 KB)
1800-175_17_ws.root (9.7 KB)
InterpolateShapes.py (5.1 KB)

Thanks a lot!

Hi @ammitra!

Looking at your code, it seems to me very likely that the problems are:

  1. You reuse the same names for different variables (some RooFit features don’t play well with this)
  2. You have different instances of the same mathematical variable, resulting in inconsistent computation graphs

The problem comes from using these external workspaces to get your pdfs:

    f1 = ROOT.TFile.Open(f'./{mT}-{m1}_{year}_ws.root')
    f2 = ROOT.TFile.Open(f'./{mT}-{m2}_{year}_ws.root')
    w1 = f1.Get('w')
    w2 = f2.Get('w')
    pdf1 = w1.pdf('model')
    pdf2 = w2.pdf('model')

    pdf1.Print("t")
    pdf2.Print("t")

The pdf objects pdf1 and pdf2, as well as their internal variables like muX, wdX etc. are actually different objects but share the same RooAbsArg names. This is discouraged, and probably there is some deduplication going on in the RooMomentMorphFuncND that makes it think that both templates should be identical and that’s why the interpolation looks like one of them.

Then, about my second point: you have different instances of mPhi and mT. The RooRealVars that you instantiate in the script, and then the objects that already exists in the workspaces in your input files. So the computation graph is inconsistent, because RooFit doesn’t know that these variables are actually the same (at least not in the general case. In practice, there is deduplication in plotting and likelihood creation that deals with this).

Can you try to suffix the names of the input models, and also their parameters so they are unique (but not the observables mPhi and mT, they are common!), and see if it works?

So that the models look like this maybe (here I suffixed with the mass):

0x2d26a530 RooProdPdf::model_125 = 7.58494e-06 [Auto,Dirty]
  0x2d208170/V- RooCrystalBall::DSCBX__copy_125 = 0.00113231 [Auto,Dirty]
    0x2cca6120/V- RooRealVar::mPhi = 310
    0x2cc306b0/V- RooRealVar::muX_125 = 122.147 +/- 2.12946
    0x2cb9adb0/V- RooRealVar::wdX_125 = 13.9833 +/- 1.65009
    0x2c9ea7b0/V- RooRealVar::a1X_125 = 1.98001 +/- 0.396451
    0x2c943260/V- RooRealVar::p1X_125 = 1.00026 +/- 70.4732
    0x25bb0130/V- RooRealVar::a2X_125 = 2.52231 +/- 0.810184
    0x25ad3280/V- RooRealVar::p2X_125 = 1.10771 +/- 3.77807
  0x28a57590/V- RooCrystalBall::DSCBY__copy_125 = 0.00669866 [Auto,Dirty]
    0x28a58e40/V- RooRealVar::mT = 2150
    0x28a5bba0/V- RooRealVar::muY_125 = 1790.98 +/- 20.4601
    0x28a5c510/V- RooRealVar::wdY_125 = 100.233 +/- 16.9621
    0x28a5cd20/V- RooRealVar::a1Y_125 = 1 +/- 0.395276
    0x28a5d530/V- RooRealVar::p1Y_125 = 13.5612 +/- 54.5128
    0x28a5dd40/V- RooRealVar::a2Y_125 = 2.60775 +/- 1.20201
    0x28a5e550/V- RooRealVar::p2Y_125 = 1.87419 +/- 4.71457

And to make sure that everything is deduplicated in the end, you can do the same old trick of importing the final model into a new RooWorkspace and then use it from there. There is also deduplication happening when importing in a RooWorkspace.

Let me know if sanitizing the computation graph like this improves the situation!

Cheers,
Jonas

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