Analytical integrals for many observables in RooClassFactory::makePdf

Hi,

This may be a very naive question, but I’m slightly confused about the syntax to add optional expressions for analytical integrals in the RooClassFactory::makePdf method. The problem is that I’m interest in multi dimensional PDF’s, which forms are basically a product of these observables. For instance, let me take a very simple example of “ypow(x-B,2)", which is somehow similar to the tutorial rf104_classfactory.C. The documentation states that the syntax follows the idea of “observable1:integral over observable1;observable2:integral over observable2”, where the “integral over observableK” is already evaluated at the min/max values of the integration. My confusion is related to the expression provided with the integration in only 1-dimension, for instance, the integral over “y” would be "1/2pow(x-B,2)*pow(y,2)” with y evaluated at its min and max values. Since this expression has still a clear dependence on “x”, I’m wondering if this is indeed the ideal expression to provide or actually it would be better to provide a single expression integrating over both x and y. Thanks in advance for any comment.

Cheers,
Rafael

What you are doing is creating a function of x isn’t it? For all x in the domain the PDF takes on some value PDF(x)=1/2*pow(x-B,2)*pow(y,2)|_{D}. ymin<D<ymax No?

Hi,

Thanks a lot for the quick reply. That’s indeed my question. I’m fine to provide several expressions integrated only in the observables of interest. However, I was just wondering why wouldn’t be better to provide it integrated over all observables at once. To be honest, the reason why I’m “manually” setting these integrals is because when I generate some new datasets from my PDF, it takes a very long to produce them. Since the function is four dimensional and it is rather complicate, my guess is that this is caused by the algorithm trying to find empirically the maximum of the function, which requires a large number of samples. Hence, my understanding is that by providing the analytical integral of the function this should be significantly faster to generate, and naively, I would expect to provide the integral already evaluated in all the observables and not individually. Am I’m missing something obvious here?

Thanks again!
Rafael

I can imagine it being quite common to want to represent the PDF along 1 dimension, e.g. x.

Regarding your second to last statement. If you want to find the maximum of the function you can’t provide an integral already evaluated over all dimensions, cause that’s just a real or complex valued number, which I am sure you are aware of. So I must be misunderstanding you?

Now, mathematically speaking, If you want to find extreme values of a four dimensional function analytically I guess you’ll have to look for stationary points etc.

I don’t understand the problem here.

Hi,

Sorry for the mis-understanding. Probably it is easier to stick to the issue that raised my question. I’ve attached one very simple script that creates a PDF dependent on 4 observables in which I provide explicitly the analytical integral as suggested. This is just a snapshot of my code where I manually set the input values. The issue that when I run “python scriptROOTForum_integral.py” I get a segmentation violation complaining:

In file included from input_line_12:9:
/home/hep/rsilvaco/Analysis/Studies/PIDCalibMoM/RooCombLegPdf.cxx:72:160: warning: expression result unused [-Wunused-value]
   if (code==1) { return (-0.0212930013306*t*y*z-0.062057257126*t*y+0.0040561212905*t*z+0.0362462619166*t-0.0256422523986*y*z-0.053579463655*y+0.056583464482*z+0.125,y:-0.01765050287334*t*x*z-0.084...
                          ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~^~~~~~
/home/hep/rsilvaco/Analysis/Studies/PIDCalibMoM/RooCombLegPdf.cxx:72:168: error: expected ')'
   if (code==1) { return (-0.0212930013306*t*y*z-0.062057257126*t*y+0.0040561212905*t*z+0.0362462619166*t-0.0256422523986*y*z-0.053579463655*y+0.056583464482*z+0.125,y:-0.01765050287334*t*x*z-0.084...
                                                                                                                                                                       ^
/home/hep/rsilvaco/Analysis/Studies/PIDCalibMoM/RooCombLegPdf.cxx:72:26: note: to match this '('
   if (code==1) { return (-0.0212930013306*t*y*z-0.062057257126*t*y+0.0040561212905*t*z+0.0362462619166*t-0.0256422523986*y*z-0.053579463655*y+0.056583464482*z+0.125,y:-0.01765050287334*t*x*z-0.084...
                         ^
Error in <ACLiC>: Dictionary generation failed!

 *** Break *** segmentation violation
 
 ===========================================================
There was a crash.
This is the entire stack trace of all threads:
===========================================================

Thread 2 (Thread 0x7f5a84423700 (LWP 9988)):
#0  0x00007f5a99878ff0 in sem_wait () from /lib64/libpthread.so.0
#1  0x00007f5a99b5da89 in PyThread_acquire_lock () from /usr/lib64/libpython2.7.so.1.0
#2  0x00007f5a99b83410 in PyEval_RestoreThread () from /usr/lib64/libpython2.7.so.1.0
#3  0x00007f5a9680a21a in ?? () from /usr/lib64/python2.7/lib-dynload/time.so
#4  0x00007f5a99b4ed75 in PyEval_EvalFrameEx () from /usr/lib64/libpython2.7.so.1.0
#5  0x00007f5a99b54961 in PyEval_EvalCodeEx () from /usr/lib64/libpython2.7.so.1.0
#6  0x00007f5a99b38eab in ?? () from /usr/lib64/libpython2.7.so.1.0
#7  0x00007f5a99b34086 in PyObject_Call () from /usr/lib64/libpython2.7.so.1.0
#8  0x00007f5a99b4fcd8 in PyEval_EvalFrameEx () from /usr/lib64/libpython2.7.so.1.0
#9  0x00007f5a99b50b6a in PyEval_EvalFrameEx () from /usr/lib64/libpython2.7.so.1.0
#10 0x00007f5a99b50b6a in PyEval_EvalFrameEx () from /usr/lib64/libpython2.7.so.1.0
#11 0x00007f5a99b54961 in PyEval_EvalCodeEx () from /usr/lib64/libpython2.7.so.1.0
#12 0x00007f5a99b38c2f in ?? () from /usr/lib64/libpython2.7.so.1.0
#13 0x00007f5a99b34086 in PyObject_Call () from /usr/lib64/libpython2.7.so.1.0
#14 0x00007f5a99b34f7a in ?? () from /usr/lib64/libpython2.7.so.1.0
#15 0x00007f5a99b34086 in PyObject_Call () from /usr/lib64/libpython2.7.so.1.0
#16 0x00007f5a99b4d630 in PyEval_CallObjectWithKeywords () from /usr/lib64/libpython2.7.so.1.0
#17 0x00007f5a99b25be8 in ?? () from /usr/lib64/libpython2.7.so.1.0
#18 0x00007f5a998730a4 in start_thread () from /lib64/libpthread.so.0
#19 0x00007f5a995a7fed in clone () from /lib64/libc.so.6

Thread 1 (Thread 0x7f5a9a009700 (LWP 9968)):
#0  0x00007f5a995799a9 in waitpid () from /lib64/libc.so.6
#1  0x00007f5a9950411b in do_system () from /lib64/libc.so.6
#2  0x00007f5a90403769 in TUnixSystem::StackTrace() () from /app/cern/root_v6.06.00/lib/libCore.so.6.06
#3  0x00007f5a90405dec in TUnixSystem::DispatchSignals(ESignals) () from /app/cern/root_v6.06.00/lib/libCore.so.6.06
#4  <signal handler called>
#5  0x00007f5a9987a73b in raise () from /lib64/libpthread.so.0
#6  0x00007f5a81e2fe78 in RooClassFactory::makePdfInstance(char const*, char const*, char const*, RooArgList const&, char const*) () from /app/cern/root_v6.06.00/lib/libRooFitCore.so.6.06
#7  0x00007f5a81e30072 in RooClassFactory::makePdfInstance(char const*, char const*, RooArgList const&, char const*) () from /app/cern/root_v6.06.00/lib/libRooFitCore.so.6.06
#8  0x00007f5a99f70097 in ?? ()
#9  0x00007ffcd5b9ad00 in ?? ()
#10 0x00007ffcd5b9a918 in ?? ()
#11 0x00007ffcd5b9a850 in ?? ()
#12 0x0000000400000000 in ?? ()
#13 0x0000000000000000 in ?? ()
===========================================================


The lines below might hint at the cause of the crash.
If they do not help you then please submit a bug report at
http://root.cern.ch/bugs. Please post the ENTIRE stack trace
from above as an attachment in addition to anything else
that might help us fixing this issue.
===========================================================
#5  0x00007f5a9987a73b in raise () from /lib64/libpthread.so.0
#6  0x00007f5a81e2fe78 in RooClassFactory::makePdfInstance(char const*, char const*, char const*, RooArgList const&, char const*) () from /app/cern/root_v6.06.00/lib/libRooFitCore.so.6.06
#7  0x00007f5a81e30072 in RooClassFactory::makePdfInstance(char const*, char const*, RooArgList const&, char const*) () from /app/cern/root_v6.06.00/lib/libRooFitCore.so.6.06
#8  0x00007f5a99f70097 in ?? ()
#9  0x00007ffcd5b9ad00 in ?? ()
#10 0x00007ffcd5b9a918 in ?? ()
#11 0x00007ffcd5b9a850 in ?? ()
#12 0x0000000400000000 in ?? ()
#13 0x0000000000000000 in ?? ()
===========================================================


Traceback (most recent call last):
  File "scriptROOTForum_integral.py", line 45, in <module>
    fitLegPol = ROOT.RooClassFactory.makePdfInstance("CombLeg", str(_simplifiedExpression), RooArgList(x,y,z,t), analyticalIntegral) 
TypeError: none of the 2 overloaded methods succeeded. Full details:
  static RooAbsPdf* RooClassFactory::makePdfInstance(const char* className, const char* name, const char* expression, const RooArgList& vars, const char* intExpression = 0) =>
    could not convert argument 3 (expected string or Unicode object, RooArgList found)
  static RooAbsPdf* RooClassFactory::makePdfInstance(const char* name, const char* expression, const RooArgList& vars, const char* intExpression = 0) =>
    problem in C++; program state has been reset

There is anything obvious I’m missing in here?

Thanks again!
scriptROOTForum_integral.py (2.5 KB)

Sorry for being slightly annoying, but there is any further suggestion about this topic?

Thanks a lot!
Rafael

Hi all again,

Unfortunately this issue has still not been solved. For simplicity, I’ve modified my previous script to a 2-variables problem to simplify the discussion. However, the issue still remains and I see a segmentation fault:

Info in <TUnixSystem::ACLiC>: creating shared library /home/hep/rsilvaco/Analysis/Studies/PIDCalibMoM/RooCombLegPdf_cxx.so
In file included from input_line_12:9:
/home/hep/rsilvaco/Analysis/Studies/PIDCalibMoM/RooCombLegPdf.cxx:66:57: warning: expression result unused [-Wunused-value]
   if (code==1) { return (0.5*pow(x.max(rangeName),2)*y - 0.5*pow(x.min(rangeName),2)*y,y:0.5*x*pow(y.max(rangeName),2) - 0.5*x*pow(y.min(rangeName),2)) ; } 
                          ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ ^ ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
/home/hep/rsilvaco/Analysis/Studies/PIDCalibMoM/RooCombLegPdf.cxx:66:90: error: expected ')'
   if (code==1) { return (0.5*pow(x.max(rangeName),2)*y - 0.5*pow(x.min(rangeName),2)*y,y:0.5*x*pow(y.max(rangeName),2) - 0.5*x*pow(y.min(rangeName),2)) ; } 
                                                                                         ^
/home/hep/rsilvaco/Analysis/Studies/PIDCalibMoM/RooCombLegPdf.cxx:66:26: note: to match this '('
   if (code==1) { return (0.5*pow(x.max(rangeName),2)*y - 0.5*pow(x.min(rangeName),2)*y,y:0.5*x*pow(y.max(rangeName),2) - 0.5*x*pow(y.min(rangeName),2)) ; } 
                         ^
Error in <ACLiC>: Dictionary generation failed!
 *** Break *** segmentation violation

Any suggestion on how to solve this? Or there is any bug in the code itself? Thanks in advance for any insight!

Cheers, Rafael
scriptROOTForum_integral_simplest.py (1.53 KB)