Logo ROOT  
Reference Guide
RooProdPdf.cxx
Go to the documentation of this file.
1 /*****************************************************************************
2  * Project: RooFit *
3  * Package: RooFitCore *
4  * @(#)root/roofitcore:$Id$
5  * Authors: *
6  * WV, Wouter Verkerke, UC Santa Barbara, verkerke@slac.stanford.edu *
7  * DK, David Kirkby, UC Irvine, dkirkby@uci.edu *
8  * *
9  * Copyright (c) 2000-2005, Regents of the University of California *
10  * and Stanford University. All rights reserved. *
11  * *
12  * Redistribution and use in source and binary forms, *
13  * with or without modification, are permitted according to the terms *
14  * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
15  *****************************************************************************/
16 
17 /**
18 \file RooProdPdf.cxx
19 \class RooProdPdf
20 \ingroup Roofitcore
21 
22 RooProdPdf is an efficient implementation of a product of PDFs of the form
23 \f[ \prod_{i=1}^{N} \mathrm{PDF}_i (x, \ldots) \f]
24 
25 PDFs may share observables. If that is the case any irreducible subset
26 of PDFs that share observables will be normalised with explicit numeric
27 integration as any built-in normalisation will no longer be valid.
28 
29 Alternatively, products using conditional PDFs can be defined, *e.g.*
30 
31 \f[ F(x|y) \cdot G(y), \f]
32 
33 meaning a PDF \f$ F(x) \f$ **given** \f$ y \f$ and a PDF \f$ G(y) \f$.
34 In this construction, \f$ F \f$ is only
35 normalised w.r.t \f$ x\f$, and \f$ G \f$ is normalised w.r.t \f$ y \f$. The product in this construction
36 is properly normalised.
37 
38 If exactly one of the component PDFs supports extended likelihood fits, the
39 product will also be usable in extended mode, returning the number of expected
40 events from the extendable component PDF. The extendable component does not
41 have to appear in any specific place in the list.
42 **/
43 
44 #include "RooProdPdf.h"
45 #include "RooRealProxy.h"
46 #include "RooProdGenContext.h"
47 #include "RooGenProdProj.h"
48 #include "RooProduct.h"
49 #include "RooNameReg.h"
50 #include "RooMsgService.h"
51 #include "RooFormulaVar.h"
52 #include "RooRealVar.h"
53 #include "RooAddition.h"
54 #include "RooGlobalFunc.h"
55 #include "RooConstVar.h"
56 #include "RooWorkspace.h"
57 #include "RooRangeBoolean.h"
58 #include "RooCustomizer.h"
59 #include "RooRealIntegral.h"
60 #include "RooTrace.h"
61 #include "RooBatchCompute.h"
62 #include "RooHelpers.h"
63 #include "strtok.h"
64 
65 #include <cstring>
66 #include <sstream>
67 #include <algorithm>
68 
69 #ifndef _WIN32
70 #include <strings.h>
71 #endif
72 
73 using namespace std;
74 
76 
77 
78 ////////////////////////////////////////////////////////////////////////////////
79 /// Default constructor
80 
82  _cutOff(0),
83  _extendedIndex(-1),
84  _useDefaultGen(kFALSE),
85  _refRangeName(0),
86  _selfNorm(kTRUE)
87 {
88  // Default constructor
90 }
91 
92 
93 
94 ////////////////////////////////////////////////////////////////////////////////
95 /// Dummy constructor
96 
97 RooProdPdf::RooProdPdf(const char *name, const char *title, Double_t cutOff) :
98  RooAbsPdf(name,title),
99  _cacheMgr(this,10),
100  _genCode(10),
101  _cutOff(cutOff),
102  _pdfList("!pdfs","List of PDFs",this),
103  _extendedIndex(-1),
104  _useDefaultGen(kFALSE),
105  _refRangeName(0),
106  _selfNorm(kTRUE)
107 {
109 }
110 
111 
112 
113 ////////////////////////////////////////////////////////////////////////////////
114 /// Constructor with 2 PDFs (most frequent use case).
115 ///
116 /// The optional cutOff parameter can be used as a speed optimization if
117 /// one or more of the PDF have sizable regions with very small values,
118 /// which would pull the entire product of PDFs to zero in those regions.
119 ///
120 /// After each PDF multiplication, the running product is compared with
121 /// the cutOff parameter. If the running product is smaller than the
122 /// cutOff value, the product series is terminated and remaining PDFs
123 /// are not evaluated.
124 ///
125 /// There is no magic value of the cutOff, the user should experiment
126 /// to find the appropriate balance between speed and precision.
127 /// If a cutoff is specified, the PDFs most likely to be small should
128 /// be put first in the product. The default cutOff value is zero.
129 ///
130 
131 RooProdPdf::RooProdPdf(const char *name, const char *title,
132  RooAbsPdf& pdf1, RooAbsPdf& pdf2, Double_t cutOff) :
133  RooAbsPdf(name,title),
134  _cacheMgr(this,10),
135  _genCode(10),
136  _cutOff(cutOff),
137  _pdfList("!pdfs","List of PDFs",this),
138  _extendedIndex(-1),
139  _useDefaultGen(kFALSE),
140  _refRangeName(0),
141  _selfNorm(kTRUE)
142 {
143  _pdfList.add(pdf1) ;
144  RooArgSet* nset1 = new RooArgSet("nset") ;
145  _pdfNSetList.Add(nset1) ;
146  if (pdf1.canBeExtended()) {
147  _extendedIndex = _pdfList.index(&pdf1) ;
148  }
149 
150  _pdfList.add(pdf2) ;
151  RooArgSet* nset2 = new RooArgSet("nset") ;
152  _pdfNSetList.Add(nset2) ;
153 
154  if (pdf2.canBeExtended()) {
155  if (_extendedIndex>=0) {
156  // Protect against multiple extended terms
157  coutW(InputArguments) << "RooProdPdf::RooProdPdf(" << GetName()
158  << ") multiple components with extended terms detected,"
159  << " product will not be extendible." << endl ;
160  _extendedIndex=-1 ;
161  } else {
162  _extendedIndex=_pdfList.index(&pdf2) ;
163  }
164  }
166 }
167 
168 
169 
170 ////////////////////////////////////////////////////////////////////////////////
171 /// Constructor from a list of PDFs.
172 ///
173 /// The optional cutOff parameter can be used as a speed optimization if
174 /// one or more of the PDF have sizable regions with very small values,
175 /// which would pull the entire product of PDFs to zero in those regions.
176 ///
177 /// After each PDF multiplication, the running product is compared with
178 /// the cutOff parameter. If the running product is smaller than the
179 /// cutOff value, the product series is terminated and remaining PDFs
180 /// are not evaluated.
181 ///
182 /// There is no magic value of the cutOff, the user should experiment
183 /// to find the appropriate balance between speed and precision.
184 /// If a cutoff is specified, the PDFs most likely to be small should
185 /// be put first in the product. The default cutOff value is zero.
186 
187 RooProdPdf::RooProdPdf(const char* name, const char* title, const RooArgList& inPdfList, Double_t cutOff) :
188  RooAbsPdf(name,title),
189  _cacheMgr(this,10),
190  _genCode(10),
191  _cutOff(cutOff),
192  _pdfList("!pdfs","List of PDFs",this),
193  _extendedIndex(-1),
194  _useDefaultGen(kFALSE),
195  _refRangeName(0),
196  _selfNorm(kTRUE)
197 {
198  RooFIter iter = inPdfList.fwdIterator();
199  RooAbsArg* arg ;
200  Int_t numExtended(0) ;
201  while((arg=(RooAbsArg*)iter.next())) {
202  RooAbsPdf* pdf = dynamic_cast<RooAbsPdf*>(arg) ;
203  if (!pdf) {
204  coutW(InputArguments) << "RooProdPdf::RooProdPdf(" << GetName() << ") list arg "
205  << arg->GetName() << " is not a PDF, ignored" << endl ;
206  continue ;
207  }
208  _pdfList.add(*pdf) ;
209 
210  RooArgSet* nset = new RooArgSet("nset") ;
211  _pdfNSetList.Add(nset) ;
212 
213  if (pdf->canBeExtended()) {
214  _extendedIndex = _pdfList.index(pdf) ;
215  numExtended++ ;
216  }
217  }
218 
219  // Protect against multiple extended terms
220  if (numExtended>1) {
221  coutW(InputArguments) << "RooProdPdf::RooProdPdf(" << GetName()
222  << ") WARNING: multiple components with extended terms detected,"
223  << " product will not be extendible." << endl ;
224  _extendedIndex = -1 ;
225  }
226 
228 }
229 
230 
231 
232 ////////////////////////////////////////////////////////////////////////////////
233 /// Constructor from named argument list.
234 /// \param[in] name Name used by RooFit
235 /// \param[in] title Title used for plotting
236 /// \param[in] fullPdfSet Set of "regular" PDFs that are normalised over all their observables
237 /// \param[in] argX Optional arguments according to table below.
238 ///
239 /// <table>
240 /// <tr><th> Argument <th> Description
241 /// <tr><td> `Conditional(pdfSet,depSet,depsAreCond=false)` <td> Add PDF to product with condition that it
242 /// only be normalized over specified observables. Any remaining observables will be conditional observables.
243 /// (Setting `depsAreCond` to true inverts this, so the observables in depSet will be the conditional observables.)
244 /// </table>
245 ///
246 /// For example, given a PDF \f$ F(x,y) \f$ and \f$ G(y) \f$,
247 ///
248 /// `RooProdPdf("P", "P", G, Conditional(F,x))` will construct a 2-dimensional PDF as follows:
249 /// \f[
250 /// P(x,y) = \frac{G(y)}{\int_y G(y)} \cdot \frac{F(x,y)}{\int_x F(x,y)},
251 /// \f]
252 ///
253 /// which is a well normalised and properly defined PDF, but different from
254 /// \f[
255 /// P'(x,y) = \frac{F(x,y) \cdot G(y)}{\int_x\int_y F(x,y) \cdot G(y)}.
256 /// \f]
257 ///
258 /// In the former case, the \f$ y \f$ distribution of \f$ P \f$ is identical to that of \f$ G \f$, while
259 /// \f$ F \f$ only is used to determine the correlation between \f$ X \f$ and \f$ Y \f$. In the latter
260 /// case, the \f$ Y \f$ distribution is defined by the product of \f$ F \f$ and \f$ G \f$.
261 ///
262 /// This \f$ P(x,y) \f$ construction is analoguous to generating events from \f$ F(x,y) \f$ with
263 /// a prototype dataset sampled from \f$ G(y) \f$.
264 
265 RooProdPdf::RooProdPdf(const char* name, const char* title, const RooArgSet& fullPdfSet,
266  const RooCmdArg& arg1, const RooCmdArg& arg2,
267  const RooCmdArg& arg3, const RooCmdArg& arg4,
268  const RooCmdArg& arg5, const RooCmdArg& arg6,
269  const RooCmdArg& arg7, const RooCmdArg& arg8) :
270  RooAbsPdf(name,title),
271  _cacheMgr(this,10),
272  _genCode(10),
273  _cutOff(0),
274  _pdfList("!pdfs","List of PDFs",this),
275  _extendedIndex(-1),
276  _useDefaultGen(kFALSE),
277  _refRangeName(0),
278  _selfNorm(kTRUE)
279 {
280  RooLinkedList l ;
281  l.Add((TObject*)&arg1) ; l.Add((TObject*)&arg2) ;
282  l.Add((TObject*)&arg3) ; l.Add((TObject*)&arg4) ;
283  l.Add((TObject*)&arg5) ; l.Add((TObject*)&arg6) ;
284  l.Add((TObject*)&arg7) ; l.Add((TObject*)&arg8) ;
285 
286  initializeFromCmdArgList(fullPdfSet,l) ;
288 }
289 
290 
291 
292 ////////////////////////////////////////////////////////////////////////////////
293 /// Constructor from named argument list
294 
295 RooProdPdf::RooProdPdf(const char* name, const char* title,
296  const RooCmdArg& arg1, const RooCmdArg& arg2,
297  const RooCmdArg& arg3, const RooCmdArg& arg4,
298  const RooCmdArg& arg5, const RooCmdArg& arg6,
299  const RooCmdArg& arg7, const RooCmdArg& arg8) :
300  RooAbsPdf(name,title),
301  _cacheMgr(this,10),
302  _genCode(10),
303  _cutOff(0),
304  _pdfList("!pdfList","List of PDFs",this),
305  _extendedIndex(-1),
306  _useDefaultGen(kFALSE),
307  _refRangeName(0),
308  _selfNorm(kTRUE)
309 {
310  RooLinkedList l ;
311  l.Add((TObject*)&arg1) ; l.Add((TObject*)&arg2) ;
312  l.Add((TObject*)&arg3) ; l.Add((TObject*)&arg4) ;
313  l.Add((TObject*)&arg5) ; l.Add((TObject*)&arg6) ;
314  l.Add((TObject*)&arg7) ; l.Add((TObject*)&arg8) ;
315 
318 }
319 
320 
321 
322 ////////////////////////////////////////////////////////////////////////////////
323 /// Internal constructor from list of named arguments
324 
325 RooProdPdf::RooProdPdf(const char* name, const char* title, const RooArgSet& fullPdfSet, const RooLinkedList& cmdArgList) :
326  RooAbsPdf(name,title),
327  _cacheMgr(this,10),
328  _genCode(10),
329  _cutOff(0),
330  _pdfList("!pdfs","List of PDFs",this),
331  _extendedIndex(-1),
332  _useDefaultGen(kFALSE),
333  _refRangeName(0),
334  _selfNorm(kTRUE)
335 {
336  initializeFromCmdArgList(fullPdfSet, cmdArgList) ;
338 }
339 
340 
341 
342 ////////////////////////////////////////////////////////////////////////////////
343 /// Copy constructor
344 
345 RooProdPdf::RooProdPdf(const RooProdPdf& other, const char* name) :
346  RooAbsPdf(other,name),
347  _cacheMgr(other._cacheMgr,this),
348  _genCode(other._genCode),
349  _cutOff(other._cutOff),
350  _pdfList("!pdfs",this,other._pdfList),
351  _extendedIndex(other._extendedIndex),
352  _useDefaultGen(other._useDefaultGen),
353  _refRangeName(other._refRangeName),
354  _selfNorm(other._selfNorm),
355  _defNormSet(other._defNormSet)
356 {
357  // Clone contents of normalizarion set list
358  RooFIter iter = other._pdfNSetList.fwdIterator();
359  RooArgSet* nset ;
360  while((nset=(RooArgSet*)iter.next())) {
361  RooArgSet* tmp = (RooArgSet*) nset->snapshot() ;
362  tmp->setName(nset->GetName()) ;
363  _pdfNSetList.Add(tmp) ;
364  }
366 }
367 
368 
369 
370 ////////////////////////////////////////////////////////////////////////////////
371 /// Initialize RooProdPdf configuration from given list of RooCmdArg configuration arguments
372 /// and set of 'regular' p.d.f.s in product
373 
375 {
376  Int_t numExtended(0) ;
377 
378  // Process set of full PDFS
379  RooFIter siter = fullPdfSet.fwdIterator() ;
380  RooAbsPdf* pdf ;
381  while((pdf=(RooAbsPdf*)siter.next())) {
382  _pdfList.add(*pdf) ;
383  RooArgSet* nset1 = new RooArgSet("nset") ;
384  _pdfNSetList.Add(nset1) ;
385 
386  if (pdf->canBeExtended()) {
387  _extendedIndex = _pdfList.index(pdf) ;
388  numExtended++ ;
389  }
390 
391  }
392 
393  // Process list of conditional PDFs
394  RooFIter iter = l.fwdIterator();
395  RooCmdArg* carg ;
396  while((carg=(RooCmdArg*)iter.next())) {
397 
398  if (0 == strcmp(carg->GetName(), "Conditional")) {
399 
400  Int_t argType = carg->getInt(0) ;
401  RooArgSet* pdfSet = (RooArgSet*) carg->getSet(0) ;
402  RooArgSet* normSet = (RooArgSet*) carg->getSet(1) ;
403 
404  RooFIter siter2 = pdfSet->fwdIterator() ;
405  RooAbsPdf* thePdf ;
406  while ((thePdf=(RooAbsPdf*)siter2.next())) {
407  _pdfList.add(*thePdf) ;
408 
409  RooArgSet* tmp = (RooArgSet*) normSet->snapshot() ;
410  tmp->setName(0 == argType ? "nset" : "cset") ;
411  _pdfNSetList.Add(tmp) ;
412 
413  if (thePdf->canBeExtended()) {
414  _extendedIndex = _pdfList.index(thePdf) ;
415  numExtended++ ;
416  }
417 
418  }
419 
420  } else if (0 != strlen(carg->GetName())) {
421  coutW(InputArguments) << "Unknown arg: " << carg->GetName() << endl ;
422  }
423  }
424 
425  // Protect against multiple extended terms
426  if (numExtended>1) {
427  coutW(InputArguments) << "RooProdPdf::RooProdPdf(" << GetName()
428  << ") WARNING: multiple components with extended terms detected,"
429  << " product will not be extendible." << endl ;
430  _extendedIndex = -1 ;
431  }
432 
433 
434 }
435 
436 
437 
438 ////////////////////////////////////////////////////////////////////////////////
439 /// Destructor
440 
442 {
443  _pdfNSetList.Delete() ;
445 }
446 
447 
449  int code ;
450  auto cache = static_cast<CacheElem*>(_cacheMgr.getObj(nset, 0, &code)) ;
451 
452  // If cache doesn't have our configuration, recalculate here
453  if (!cache) {
454  code = getPartIntList(nset, nullptr) ;
455  cache = static_cast<CacheElem*>(_cacheMgr.getObj(nset, 0, &code)) ;
456  }
457  return cache;
458 }
459 
460 
461 ////////////////////////////////////////////////////////////////////////////////
462 /// Calculate current value of object
463 
465 {
466  return calculate(*getCacheElem(_normSet)) ;
467 }
468 
469 
470 
471 ////////////////////////////////////////////////////////////////////////////////
472 /// Calculate running product of pdfs terms, using the supplied
473 /// normalization set in 'normSetList' for each component
474 
476 {
477  //cout << "RooProdPdf::calculate from cache" << endl ;
478 
479  if (cache._isRearranged) {
480  if (dologD(Eval)) {
481  cxcoutD(Eval) << "RooProdPdf::calculate(" << GetName() << ") rearranged product calculation"
482  << " calculate: num = " << cache._rearrangedNum->GetName() << " = " << cache._rearrangedNum->getVal() << endl ;
483 // cache._rearrangedNum->printComponentTree("",0,5) ;
484  cxcoutD(Eval) << "calculate: den = " << cache._rearrangedDen->GetName() << " = " << cache._rearrangedDen->getVal() << endl ;
485 // cache._rearrangedDen->printComponentTree("",0,5) ;
486  }
487 
488  return cache._rearrangedNum->getVal() / cache._rearrangedDen->getVal();
489  } else {
490 
491  Double_t value = 1.0;
492  assert(cache._normList.size() == cache._partList.size());
493  for (std::size_t i = 0; i < cache._partList.size(); ++i) {
494  const auto& partInt = static_cast<const RooAbsReal&>(cache._partList[i]);
495  const auto normSet = cache._normList[i].get();
496 
497  const Double_t piVal = partInt.getVal(normSet->getSize() > 0 ? normSet : nullptr);
498  value *= piVal ;
499  if (value <= _cutOff) break;
500  }
501 
502  return value ;
503  }
504 }
505 
506 
507 ////////////////////////////////////////////////////////////////////////////////
508 /// Evaluate product of PDFs using input data in `evalData`.
510  int code;
511  auto cache = static_cast<CacheElem*>(_cacheMgr.getObj(normSet, nullptr, &code));
512 
513  // If cache doesn't have our configuration, recalculate here
514  if (!cache) {
515  code = getPartIntList(normSet, nullptr);
516  cache = static_cast<CacheElem*>(_cacheMgr.getObj(normSet, nullptr, &code));
517  }
518 
519  if (cache->_isRearranged) {
520  auto numerator = cache->_rearrangedNum->getValues(evalData, normSet);
521  auto denominator = cache->_rearrangedDen->getValues(evalData, normSet);
522  auto outputs = evalData.makeBatch(this, numerator.size());
523 
524  for (std::size_t i=0; i < outputs.size(); ++i) {
525  outputs[i] = numerator[i] / denominator[i];
526  }
527 
528  return outputs;
529  } else {
530  assert(cache->_normList.size() == cache->_partList.size());
531  RooSpan<double> outputs;
532  for (std::size_t i = 0; i < cache->_partList.size(); ++i) {
533  const auto& partInt = static_cast<const RooAbsReal&>(cache->_partList[i]);
534  const auto partNorm = cache->_normList[i].get();
535 
536  const auto partialInt = partInt.getValues(evalData, partNorm->getSize() > 0 ? partNorm : nullptr);
537 
538  if (outputs.empty()) {
539  outputs = evalData.makeBatch(this, partialInt.size());
540  std::fill(outputs.begin(), outputs.end(), 1.);
541  } else if (outputs.size() == 1 && partialInt.size() > 1) {
542  const double val = outputs[0];
543  outputs = evalData.makeBatch(this, partialInt.size());
544  std::fill(outputs.begin(), outputs.end(), val);
545  }
546  assert(outputs.size() == partialInt.size());
547 
548  for (std::size_t j=0; j < outputs.size(); ++j) {
549  outputs[j] *= partialInt[j];
550  }
551  }
552 
553  return outputs;
554  }
555 }
556 
557 namespace {
558 
559 template<class T>
560 void eraseNullptrs(std::vector<T*>& v) {
561  v.erase(std::remove_if(v.begin(), v.end(), [](T* x){ return x == nullptr; } ), v.end());
562 }
563 
564 void removeCommon(std::vector<RooAbsArg*> &v, std::span<RooAbsArg * const> other) {
565 
566  for (auto const& arg : other) {
567  auto namePtrMatch = [&arg](const RooAbsArg* elm) {
568  return elm->namePtr() == arg->namePtr();
569  };
570 
571  auto found = std::find_if(v.begin(), v.end(), namePtrMatch);
572  if(found != v.end()) {
573  *found = nullptr;
574  }
575  }
576  eraseNullptrs(v);
577 }
578 
579 void addCommon(std::vector<RooAbsArg*> &v, std::vector<RooAbsArg*> const& o1, std::vector<RooAbsArg*> const& o2) {
580 
581  for (auto const& arg : o1) {
582  auto namePtrMatch = [&arg](const RooAbsArg* elm) {
583  return elm->namePtr() == arg->namePtr();
584  };
585 
586  if(std::find_if(o2.begin(), o2.end(), namePtrMatch) != o2.end()) {
587  v.push_back(arg);
588  }
589  }
590 }
591 
592 }
593 
594 
595 ////////////////////////////////////////////////////////////////////////////////
596 /// Factorize product in irreducible terms for given choice of integration/normalization
597 
598 void RooProdPdf::factorizeProduct(const RooArgSet& normSet, const RooArgSet& intSet,
599  RooLinkedList& termList, RooLinkedList& normList,
600  RooLinkedList& impDepList, RooLinkedList& crossDepList,
601  RooLinkedList& intList) const
602 {
603  // List of all term dependents: normalization and imported
604  std::vector<RooArgSet> depAllList;
605  std::vector<RooArgSet> depIntNoNormList;
606 
607  // Setup lists for factorization terms and their dependents
608  RooArgSet* term(0);
609  RooArgSet* termNormDeps(0);
610  RooArgSet* termIntDeps(0);
611  RooArgSet* termIntNoNormDeps(0);
612 
613  std::vector<RooAbsArg*> pdfIntNoNormDeps;
614  std::vector<RooAbsArg*> pdfIntSet;
615  std::vector<RooAbsArg*> pdfNSet;
616  std::vector<RooAbsArg*> pdfCSet;
617  std::vector<RooAbsArg*> pdfNormDeps; // Dependents to be normalized for the PDF
618  std::vector<RooAbsArg*> pdfAllDeps; // All dependents of this PDF
619 
620  // Loop over the PDFs
621  RooAbsPdf* pdf;
622  RooArgSet* pdfNSetOrig;
623  for (RooLinkedListIter pdfIter = _pdfList.iterator(),
624  nIter = _pdfNSetList.iterator();
625  (pdfNSetOrig = static_cast<RooArgSet*>(nIter.Next()),
626  pdf = static_cast<RooAbsPdf*>(pdfIter.Next())); ) {
627  pdfNSet.clear();
628  pdfCSet.clear();
629 
630  // Make iterator over tree leaf node list to get the observables.
631  // This code is borrowed from RooAgsPdf::getObservables.
632  // RooAbsArg::treeNodeServer list is relatively expensive, so we only do it
633  // once and use it in a lambda function.
634  RooArgSet pdfLeafList("leafNodeServerList") ;
635  pdf->treeNodeServerList(&pdfLeafList,0,kFALSE,kTRUE,true) ;
636  auto getObservables = [&pdfLeafList](
637  std::vector<RooAbsArg*> & out,
638  const RooArgSet& dataList) {
639  for (const auto arg : pdfLeafList) {
640  if (arg->dependsOnValue(dataList) && arg->isLValue()) {
641  out.push_back(arg) ;
642  }
643  }
644  };
645 
646  // Reduce pdfNSet to actual dependents
647  if (0 == strcmp("nset", pdfNSetOrig->GetName())) {
648  getObservables(pdfNSet, *pdfNSetOrig);
649  } else if (0 == strcmp("cset", pdfNSetOrig->GetName())) {
650  getObservables(pdfNSet, normSet);
651  removeCommon(pdfNSet, pdfNSetOrig->get());
652  pdfCSet = pdfNSetOrig->get();
653  } else {
654  // Legacy mode. Interpret at NSet for backward compatibility
655  getObservables(pdfNSet, *pdfNSetOrig);
656  }
657 
658 
659  pdfNormDeps.clear();
660  pdfAllDeps.clear();
661 
662  // Make list of all dependents of this PDF
663  getObservables(pdfAllDeps, normSet);
664 
665 
666 // cout << GetName() << ": pdf = " << pdf->GetName() << " pdfAllDeps = " << pdfAllDeps << " pdfNSet = " << *pdfNSet << " pdfCSet = " << *pdfCSet << endl;
667 
668  // Make list of normalization dependents for this PDF;
669  if (!pdfNSet.empty()) {
670  // PDF is conditional
671  addCommon(pdfNormDeps, pdfAllDeps, pdfNSet);
672  } else {
673  // PDF is regular
674  pdfNormDeps = pdfAllDeps;
675  }
676 
677 // cout << GetName() << ": pdfNormDeps for " << pdf->GetName() << " = " << pdfNormDeps << endl;
678 
679  pdfIntSet.clear();
680  getObservables(pdfIntSet, intSet) ;
681 
682  // WVE if we have no norm deps, conditional observables should be taken out of pdfIntSet
683  if (pdfNormDeps.empty() && !pdfCSet.empty()) {
684  removeCommon(pdfIntSet, pdfCSet);
685 // cout << GetName() << ": have no norm deps, removing conditional observables from intset" << endl;
686  }
687 
688  pdfIntNoNormDeps.clear();
689  pdfIntNoNormDeps = pdfIntSet;
690  removeCommon(pdfIntNoNormDeps, pdfNormDeps);
691 
692 // cout << GetName() << ": pdf = " << pdf->GetName() << " intset = " << *pdfIntSet << " pdfIntNoNormDeps = " << pdfIntNoNormDeps << endl;
693 
694  // Check if this PDF has dependents overlapping with one of the existing terms
695  bool done = false;
696  int j = 0;
697  for (RooLinkedListIter lIter = termList.iterator(),
698  ldIter = normList.iterator();
699  (termNormDeps = static_cast<RooArgSet*>(ldIter.Next()),
700  term = static_cast<RooArgSet*>(lIter.Next())); ++j) {
701  // PDF should be added to existing term if
702  // 1) It has overlapping normalization dependents with any other PDF in existing term
703  // 2) It has overlapping dependents of any class for which integration is requested
704  // 3) If normalization happens over multiple ranges, and those ranges are both defined
705  // in either observable
706 
707  bool normOverlap = termNormDeps->overlaps(pdfNormDeps.begin(), pdfNormDeps.end());
708  //Bool_t intOverlap = pdfIntSet->overlaps(*termAllDeps);
709 
710  if (normOverlap) {
711 // cout << GetName() << ": this term overlaps with term " << (*term) << " in normalization observables" << endl;
712 
713  term->add(*pdf);
714  termNormDeps->add(pdfNormDeps.begin(), pdfNormDeps.end(), false);
715  depAllList[j].add(pdfAllDeps.begin(), pdfAllDeps.end(), false);
716  if (termIntDeps) {
717  termIntDeps->add(pdfIntSet.begin(), pdfIntSet.end(), false);
718  }
719  if (termIntNoNormDeps) {
720  termIntNoNormDeps->add(pdfIntNoNormDeps.begin(), pdfIntNoNormDeps.end(), false);
721  }
722  termIntNoNormDeps->add(pdfIntNoNormDeps.begin(), pdfIntNoNormDeps.end(), false);
723  done = true;
724  }
725  }
726 
727  // If not, create a new term
728  if (!done) {
729  if (!(pdfNormDeps.empty() && pdfAllDeps.empty() &&
730  pdfIntSet.empty()) || 0 == normSet.getSize()) {
731 // cout << GetName() << ": creating new term" << endl;
732  term = new RooArgSet("term");
733  termNormDeps = new RooArgSet("termNormDeps");
734  depAllList.emplace_back(pdfAllDeps.begin(), pdfAllDeps.end(), "termAllDeps");
735  termIntDeps = new RooArgSet(pdfIntSet.begin(), pdfIntSet.end(), "termIntDeps");
736  depIntNoNormList.emplace_back(pdfIntNoNormDeps.begin(), pdfIntNoNormDeps.end(), "termIntNoNormDeps");
737  termIntNoNormDeps = &depIntNoNormList.back();
738 
739  term->add(*pdf);
740  termNormDeps->add(pdfNormDeps.begin(), pdfNormDeps.end(), false);
741 
742  termList.Add(term);
743  normList.Add(termNormDeps);
744  intList.Add(termIntDeps);
745  }
746  }
747 
748  }
749 
750  // Loop over list of terms again to determine 'imported' observables
751  int i = 0;
752  RooArgSet *normDeps;
753  for (RooLinkedListIter lIter = termList.iterator(),
754  ldIter = normList.iterator();
755  (normDeps = static_cast<RooArgSet*>(ldIter.Next()),
756  term= static_cast<RooArgSet*>(lIter.Next())); ++i) {
757  // Make list of wholly imported dependents
758  RooArgSet impDeps(depAllList[i]);
759  impDeps.remove(*normDeps, kTRUE, kTRUE);
760  impDepList.Add(impDeps.snapshot());
761 // cout << GetName() << ": list of imported dependents for term " << (*term) << " set to " << impDeps << endl ;
762 
763  // Make list of cross dependents (term is self contained for these dependents,
764  // but components import dependents from other components)
765  auto crossDeps = std::unique_ptr<RooAbsCollection>{depIntNoNormList[i].selectCommon(*normDeps)};
766  crossDepList.Add(crossDeps->snapshot());
767 // cout << GetName() << ": list of cross dependents for term " << (*term) << " set to " << *crossDeps << endl ;
768  }
769 
770  return;
771 }
772 
773 
774 
775 
776 ////////////////////////////////////////////////////////////////////////////////
777 /// Return list of (partial) integrals of product terms for integration
778 /// of p.d.f over observables iset while normalization over observables nset.
779 /// Also return list of normalization sets to be used to evaluate
780 /// each component in the list correctly.
781 
782 Int_t RooProdPdf::getPartIntList(const RooArgSet* nset, const RooArgSet* iset, const char* isetRangeName) const
783 {
784 // cout << " FOLKERT::RooProdPdf::getPartIntList(" << GetName() <<") nset = " << (nset?*nset:RooArgSet()) << endl
785 // << " _normRange = " << _normRange << endl
786 // << " iset = " << (iset?*iset:RooArgSet()) << endl
787 // << " isetRangeName = " << (isetRangeName?isetRangeName:"<null>") << endl ;
788 
789  // Check if this configuration was created before
790  Int_t sterileIdx(-1);
791 
792  CacheElem* cache = (CacheElem*) _cacheMgr.getObj(nset,iset,&sterileIdx,isetRangeName);
793  if (cache) {
794  return _cacheMgr.lastIndex();
795  }
796 
797  // Create containers for partial integral components to be generated
798  cache = new CacheElem;
799 
800  // Factorize the product in irreducible terms for this nset
801  RooLinkedList terms, norms, imp, ints, cross;
802 // cout << "RooProdPdf::getPIL -- now calling factorizeProduct()" << endl ;
803 
804 
805  // Normalization set used for factorization
806  RooArgSet factNset(nset ? (*nset) : _defNormSet);
807 // cout << GetName() << "factNset = " << factNset << endl ;
808 
809  factorizeProduct(factNset, iset ? (*iset) : RooArgSet(), terms, norms, imp, cross, ints);
810 
811  RooArgSet *norm, *integ, *xdeps, *imps;
812 
813  // Group irriducible terms that need to be (partially) integrated together
814  std::list<std::vector<RooArgSet*>> groupedList;
815  RooArgSet outerIntDeps;
816 // cout << "RooProdPdf::getPIL -- now calling groupProductTerms()" << endl;
817  groupProductTerms(groupedList, outerIntDeps, terms, norms, imp, ints, cross);
818 
819  // Loop over groups
820 // cout<<"FK: pdf("<<GetName()<<") Starting selecting F(x|y)!"<<endl;
821  // Find groups of type F(x|y), i.e. termImpSet!=0, construct ratio object
822  map<string, RooArgSet> ratioTerms;
823  for (auto const& group : groupedList) {
824  if (1 == group.size()) {
825 // cout<<"FK: Starting Single Term"<<endl;
826 
827  RooArgSet* term = group[0];
828 
829  Int_t termIdx = terms.IndexOf(term);
830  norm=(RooArgSet*) norms.At(termIdx);
831  imps=(RooArgSet*)imp.At(termIdx);
832  RooArgSet termNSet(*norm), termImpSet(*imps);
833 
834 // cout<<"FK: termImpSet.getSize() = "<<termImpSet.getSize()<< " " << termImpSet << endl;
835 // cout<<"FK: _refRangeName = "<<_refRangeName<<endl;
836 
837  if (termImpSet.getSize() > 0 && 0 != _refRangeName) {
838 
839 // cout << "WVE now here" << endl;
840 
841  // WVE we can skip this if the ref range is equal to the normalization range
842  Bool_t rangeIdentical(kTRUE);
843  RooFIter niter = termNSet.fwdIterator();
844  RooRealVar* normObs;
845 // cout << "_normRange = " << _normRange << " _refRangeName = " << RooNameReg::str(_refRangeName) << endl ;
846  while ((normObs = (RooRealVar*) niter.next())) {
847  //FK: Here the refRange should be compared to _normRange, if it's set, and to the normObs range if it's not set
848  if (_normRange.Length() > 0) {
849  if (normObs->getMin(_normRange.Data()) != normObs->getMin(RooNameReg::str(_refRangeName))) rangeIdentical = kFALSE;
850  if (normObs->getMax(_normRange.Data()) != normObs->getMax(RooNameReg::str(_refRangeName))) rangeIdentical = kFALSE;
851  }
852  else{
853  if (normObs->getMin() != normObs->getMin(RooNameReg::str(_refRangeName))) rangeIdentical = kFALSE;
854  if (normObs->getMax() != normObs->getMax(RooNameReg::str(_refRangeName))) rangeIdentical = kFALSE;
855  }
856  }
857 // cout<<"FK: rangeIdentical Single = "<<(rangeIdentical ? 'T':'F')<<endl;
858  // coverity[CONSTANT_EXPRESSION_RESULT]
859  // LM : avoid making integral ratio if range is the same. Why was not included ??? (same at line 857)
860  if (!rangeIdentical ) {
861 // cout << "PREPARING RATIO HERE (SINGLE TERM)" << endl ;
862  RooAbsReal* ratio = makeCondPdfRatioCorr(*(RooAbsReal*)term->first(), termNSet, termImpSet, normRange(), RooNameReg::str(_refRangeName));
863  ostringstream str; termImpSet.printValue(str);
864 // cout << GetName() << "inserting ratio term" << endl;
865  ratioTerms[str.str()].add(*ratio);
866  }
867  }
868 
869  } else {
870 // cout<<"FK: Starting Composite Term"<<endl;
871 
872  RooArgSet compTermSet, compTermNorm;
873  for (auto const& term : group) {
874 
875  Int_t termIdx = terms.IndexOf(term);
876  norm=(RooArgSet*) norms.At(termIdx);
877  imps=(RooArgSet*)imp.At(termIdx);
878  RooArgSet termNSet(*norm), termImpSet(*imps);
879 
880  if (termImpSet.getSize() > 0 && 0 != _refRangeName) {
881 
882  // WVE we can skip this if the ref range is equal to the normalization range
883  Bool_t rangeIdentical(kTRUE);
884  RooFIter niter = termNSet.fwdIterator();
885  RooRealVar* normObs;
886  //FK: Here the refRange should be compared to _normRange, if it's set, and to the normObs range if it's not set
887  if(_normRange.Length() > 0) {
888  while ((normObs = (RooRealVar*) niter.next())) {
889  if (normObs->getMin(_normRange.Data()) != normObs->getMin(RooNameReg::str(_refRangeName))) rangeIdentical = kFALSE;
890  if (normObs->getMax(_normRange.Data()) != normObs->getMax(RooNameReg::str(_refRangeName))) rangeIdentical = kFALSE;
891  }
892  } else {
893  while ((normObs = (RooRealVar*) niter.next())) {
894  if (normObs->getMin() != normObs->getMin(RooNameReg::str(_refRangeName))) rangeIdentical = kFALSE;
895  if (normObs->getMax() != normObs->getMax(RooNameReg::str(_refRangeName))) rangeIdentical = kFALSE;
896  }
897  }
898 // cout<<"FK: rangeIdentical Composite = "<<(rangeIdentical ? 'T':'F') <<endl;
899  if (!rangeIdentical ) {
900 // cout << "PREPARING RATIO HERE (COMPOSITE TERM)" << endl ;
901  RooAbsReal* ratio = makeCondPdfRatioCorr(*(RooAbsReal*)term->first(), termNSet, termImpSet, normRange(), RooNameReg::str(_refRangeName));
902  ostringstream str; termImpSet.printValue(str);
903  ratioTerms[str.str()].add(*ratio);
904  }
905  }
906  }
907  }
908 
909  }
910 
911  // Find groups with y as termNSet
912  // Replace G(y) with (G(y),ratio)
913  for (auto const& group : groupedList) {
914  if (1 == group.size()) {
915  RooArgSet* term = group[0];
916 
917  Int_t termIdx = terms.IndexOf(term);
918  norm = (RooArgSet*) norms.At(termIdx);
919  imps = (RooArgSet*) imp.At(termIdx);
920  RooArgSet termNSet(*norm), termImpSet(*imps);
921 
922  // If termNset matches index of ratioTerms, insert ratio here
923  ostringstream str; termNSet.printValue(str);
924  if (ratioTerms[str.str()].getSize() > 0) {
925 // cout << "MUST INSERT RATIO OBJECT IN TERM (SINGLE) " << *term << endl;
926  term->add(ratioTerms[str.str()]);
927  }
928  } else {
929  RooArgSet compTermSet, compTermNorm;
930  for (auto const& term : group) {
931  Int_t termIdx = terms.IndexOf(term);
932  norm = (RooArgSet*) norms.At(termIdx);
933  imps = (RooArgSet*) imp.At(termIdx);
934  RooArgSet termNSet(*norm), termImpSet(*imps);
935 
936  // If termNset matches index of ratioTerms, insert ratio here
937  ostringstream str; termNSet.printValue(str);
938  if (ratioTerms[str.str()].getSize() > 0) {
939 // cout << "MUST INSERT RATIO OBJECT IN TERM (COMPOSITE)" << *term << endl;
940  term->add(ratioTerms[str.str()]);
941  }
942  }
943  }
944  }
945 
946  for (auto const& group : groupedList) {
947 // cout << GetName() << ":now processing group" << endl;
948 // group->Print("1");
949 
950  if (1 == group.size()) {
951 // cout << "processing atomic item" << endl;
952  RooArgSet* term = group[0];
953 
954  Int_t termIdx = terms.IndexOf(term);
955  norm = (RooArgSet*) norms.At(termIdx);
956  integ = (RooArgSet*) ints.At(termIdx);
957  xdeps = (RooArgSet*) cross.At(termIdx);
958  imps = (RooArgSet*) imp.At(termIdx);
959 
960  RooArgSet termNSet, termISet, termXSet, termImpSet;
961 
962  // Take list of normalization, integrated dependents from factorization algorithm
963  termISet.add(*integ);
964  termNSet.add(*norm);
965 
966  // Cross-imported integrated dependents
967  termXSet.add(*xdeps);
968  termImpSet.add(*imps);
969 
970 // cout << GetName() << ": termISet = " << termISet << endl;
971 // cout << GetName() << ": termNSet = " << termNSet << endl;
972 // cout << GetName() << ": termXSet = " << termXSet << endl;
973 // cout << GetName() << ": termImpSet = " << termImpSet << endl;
974 
975  // Add prefab term to partIntList.
976  Bool_t isOwned(kFALSE);
977  vector<RooAbsReal*> func = processProductTerm(nset, iset, isetRangeName, term, termNSet, termISet, isOwned);
978  if (func[0]) {
979  cache->_partList.add(*func[0]);
980  if (isOwned) cache->_ownedList.addOwned(*func[0]);
981 
982  cache->_normList.emplace_back(norm->snapshot(kFALSE));
983 
984  cache->_numList.addOwned(*func[1]);
985  cache->_denList.addOwned(*func[2]);
986 // cout << "func[0]=" << func[0]->IsA()->GetName() << "::" << func[0]->GetName() << endl;
987 // cout << "func[1]=" << func[1]->IsA()->GetName() << "::" << func[1]->GetName() << endl;
988 // cout << "func[2]=" << func[2]->IsA()->GetName() << "::" << func[2]->GetName() << endl;
989  }
990  } else {
991 // cout << "processing composite item" << endl;
992  RooArgSet compTermSet, compTermNorm, compTermNum, compTermDen;
993  for (auto const& term : group) {
994 // cout << GetName() << ": processing term " << (*term) << " of composite item" << endl ;
995  Int_t termIdx = terms.IndexOf(term);
996  norm = (RooArgSet*) norms.At(termIdx);
997  integ = (RooArgSet*) ints.At(termIdx);
998  xdeps = (RooArgSet*) cross.At(termIdx);
999  imps = (RooArgSet*) imp.At(termIdx);
1000 
1001  RooArgSet termNSet, termISet, termXSet, termImpSet;
1002  termISet.add(*integ);
1003  termNSet.add(*norm);
1004  termXSet.add(*xdeps);
1005  termImpSet.add(*imps);
1006 
1007  // Remove outer integration dependents from termISet
1008  termISet.remove(outerIntDeps, kTRUE, kTRUE);
1009 // cout << "termISet = "; termISet.Print("1");
1010 
1011 // cout << GetName() << ": termISet = " << termISet << endl;
1012 // cout << GetName() << ": termNSet = " << termNSet << endl;
1013 // cout << GetName() << ": termXSet = " << termXSet << endl;
1014 // cout << GetName() << ": termImpSet = " << termImpSet << endl;
1015 
1016  Bool_t isOwned = false;
1017  vector<RooAbsReal*> func = processProductTerm(nset, iset, isetRangeName, term, termNSet, termISet, isOwned, kTRUE);
1018 // cout << GetName() << ": created composite term component " << func[0]->GetName() << endl;
1019  if (func[0]) {
1020  compTermSet.add(*func[0]);
1021  if (isOwned) cache->_ownedList.addOwned(*func[0]);
1022  compTermNorm.add(*norm, kFALSE);
1023 
1024  compTermNum.add(*func[1]);
1025  compTermDen.add(*func[2]);
1026  //cache->_numList.add(*func[1]);
1027  //cache->_denList.add(*func[2]);
1028 
1029 // cout << "func[0]=" << func[0]->IsA()->GetName() << "::" << func[0]->GetName() << endl;
1030 // cout << "func[1]=" << func[1]->IsA()->GetName() << "::" << func[1]->GetName() << endl;
1031 // cout << "func[2]=" << func[2]->IsA()->GetName() << "::" << func[2]->GetName() << endl;
1032  }
1033  }
1034 
1035 // cout << GetName() << ": constructing special composite product" << endl;
1036 // cout << GetName() << ": compTermSet = " ; compTermSet.Print("1");
1037 
1038  // WVE THIS NEEDS TO BE REARRANGED
1039 
1040  // compTermset is set van partial integrals to be multiplied
1041  // prodtmp = product (compTermSet)
1042  // inttmp = int ( prodtmp ) d (outerIntDeps) _range_isetRangeName
1043 
1044  const std::string prodname = makeRGPPName("SPECPROD", compTermSet, outerIntDeps, RooArgSet(), isetRangeName);
1045  RooProduct* prodtmp = new RooProduct(prodname.c_str(), prodname.c_str(), compTermSet);
1046  cache->_ownedList.addOwned(*prodtmp);
1047 
1048  const std::string intname = makeRGPPName("SPECINT", compTermSet, outerIntDeps, RooArgSet(), isetRangeName);
1049  RooRealIntegral* inttmp = new RooRealIntegral(intname.c_str(), intname.c_str(), *prodtmp, outerIntDeps, 0, 0, isetRangeName);
1050  inttmp->setStringAttribute("PROD_TERM_TYPE", "SPECINT");
1051 
1052  cache->_ownedList.addOwned(*inttmp);
1053  cache->_partList.add(*inttmp);
1054 
1055  // Product of numerator terms
1056  const string prodname_num = makeRGPPName("SPECPROD_NUM", compTermNum, RooArgSet(), RooArgSet(), 0);
1057  RooProduct* prodtmp_num = new RooProduct(prodname_num.c_str(), prodname_num.c_str(), compTermNum);
1058  prodtmp_num->addOwnedComponents(compTermNum);
1059  cache->_ownedList.addOwned(*prodtmp_num);
1060 
1061  // Product of denominator terms
1062  const string prodname_den = makeRGPPName("SPECPROD_DEN", compTermDen, RooArgSet(), RooArgSet(), 0);
1063  RooProduct* prodtmp_den = new RooProduct(prodname_den.c_str(), prodname_den.c_str(), compTermDen);
1064  prodtmp_den->addOwnedComponents(compTermDen);
1065  cache->_ownedList.addOwned(*prodtmp_den);
1066 
1067  // Ratio
1068  string name = Form("SPEC_RATIO(%s,%s)", prodname_num.c_str(), prodname_den.c_str());
1069  RooFormulaVar* ndr = new RooFormulaVar(name.c_str(), "@0/@1", RooArgList(*prodtmp_num, *prodtmp_den));
1070 
1071  // Integral of ratio
1072  RooAbsReal* numtmp = ndr->createIntegral(outerIntDeps,isetRangeName);
1073  numtmp->addOwnedComponents(*ndr);
1074 
1075  cache->_numList.addOwned(*numtmp);
1076  cache->_denList.addOwned(*(RooAbsArg*)RooFit::RooConst(1).clone("1"));
1077  cache->_normList.emplace_back(compTermNorm.snapshot(kFALSE));
1078  }
1079  }
1080 
1081  // Store the partial integral list and return the assigned code
1082  Int_t returnCode = _cacheMgr.setObj(nset, iset, (RooAbsCacheElement*)cache, RooNameReg::ptr(isetRangeName));
1083 
1084  // WVE DEBUG PRINTING
1085 // cout << "RooProdPdf::getPartIntList(" << GetName() << ") made cache " << cache << " with the following nset pointers ";
1086 // TIterator* nliter = nsetList->MakeIterator();
1087 // RooArgSet* ns;
1088 // while((ns=(RooArgSet*)nliter->Next())) {
1089 // cout << ns << " ";
1090 // }
1091 // cout << endl;
1092 // delete nliter;
1093 
1094 // cout << " FOLKERT::RooProdPdf::getPartIntList END(" << GetName() <<") nset = " << (nset?*nset:RooArgSet()) << endl
1095 // << " _normRange = " << _normRange << endl
1096 // << " iset = " << (iset?*iset:RooArgSet()) << endl
1097 // << " partList = ";
1098 // if(partListPointer) partListPointer->Print();
1099 // cout << " nsetList = ";
1100 // if(nsetListPointer) nsetListPointer->Print("");
1101 // cout << " code = " << returnCode << endl
1102 // << " isetRangeName = " << (isetRangeName?isetRangeName:"<null>") << endl;
1103 
1104 
1105  // Need to rearrange product in case of multiple ranges
1106  if (_normRange.Contains(",")) {
1107  rearrangeProduct(*cache);
1108  }
1109 
1110  // We own contents of all lists filled by factorizeProduct()
1111  terms.Delete();
1112  ints.Delete();
1113  imp.Delete();
1114  norms.Delete();
1115  cross.Delete();
1116 
1117  return returnCode;
1118 }
1119 
1120 
1121 
1122 
1123 ////////////////////////////////////////////////////////////////////////////////
1124 /// For single normalization ranges
1125 
1126 RooAbsReal* RooProdPdf::makeCondPdfRatioCorr(RooAbsReal& pdf, const RooArgSet& termNset, const RooArgSet& /*termImpSet*/, const char* normRangeTmp, const char* refRange) const
1127 {
1128  RooAbsReal* ratio_num = pdf.createIntegral(termNset,normRangeTmp) ;
1129  RooAbsReal* ratio_den = pdf.createIntegral(termNset,refRange) ;
1130  RooFormulaVar* ratio = new RooFormulaVar(Form("ratio(%s,%s)",ratio_num->GetName(),ratio_den->GetName()),"@0/@1",
1131  RooArgList(*ratio_num,*ratio_den)) ;
1132 
1133  ratio->addOwnedComponents(RooArgSet(*ratio_num,*ratio_den)) ;
1134  ratio->setAttribute("RATIO_TERM") ;
1135  return ratio ;
1136 }
1137 
1138 
1139 
1140 
1141 ////////////////////////////////////////////////////////////////////////////////
1142 
1144 {
1145  RooAbsReal* part, *num, *den ;
1146  RooArgSet nomList ;
1147 
1148  list<string> rangeComps ;
1149  {
1150  char* buf = new char[strlen(_normRange.Data()) + 1] ;
1151  strcpy(buf,_normRange.Data()) ;
1152  char* save(0) ;
1153  char* token = R__STRTOK_R(buf,",",&save) ;
1154  while(token) {
1155  rangeComps.push_back(token) ;
1156  token = R__STRTOK_R(0,",",&save) ;
1157  }
1158  delete[] buf;
1159  }
1160 
1161 
1162  map<string,RooArgSet> denListList ;
1163  RooArgSet specIntDeps ;
1164  string specIntRange ;
1165 
1166 // cout << "THIS IS REARRANGEPRODUCT" << endl ;
1167 
1168  RooFIter iterp = cache._partList.fwdIterator() ;
1169  RooFIter iter1 = cache._numList.fwdIterator() ;
1170  RooFIter iter2 = cache._denList.fwdIterator() ;
1171  while((part=(RooAbsReal*)iterp.next())) {
1172 
1173  num = (RooAbsReal*) iter1.next() ;
1174  den = (RooAbsReal*) iter2.next() ;
1175 
1176 // cout << "now processing part " << part->GetName() << " of type " << part->getStringAttribute("PROD_TERM_TYPE") << endl ;
1177 // cout << "corresponding numerator = " << num->GetName() << endl ;
1178 // cout << "corresponding denominator = " << den->GetName() << endl ;
1179 
1180 
1181  RooFormulaVar* ratio(0) ;
1182  RooArgSet origNumTerm ;
1183 
1184  if (string("SPECINT")==part->getStringAttribute("PROD_TERM_TYPE")) {
1185 
1186  RooRealIntegral* orig = (RooRealIntegral*) num;
1187  RooFormulaVar* specratio = (RooFormulaVar*) &orig->integrand() ;
1188  RooProduct* func = (RooProduct*) specratio->getParameter(0) ;
1189 
1190  RooArgSet* comps = orig->getComponents() ;
1191  RooFIter iter = comps->fwdIterator() ;
1192  RooAbsArg* carg ;
1193  while((carg=(RooAbsArg*)iter.next())) {
1194  if (carg->getAttribute("RATIO_TERM")) {
1195  ratio = (RooFormulaVar*)carg ;
1196  break ;
1197  }
1198  }
1199  delete comps ;
1200 
1201  if (ratio) {
1202  RooCustomizer cust(*func,"blah") ;
1203  cust.replaceArg(*ratio,RooFit::RooConst(1)) ;
1204  RooAbsArg* funcCust = cust.build() ;
1205 // cout << "customized function = " << endl ;
1206 // funcCust->printComponentTree() ;
1207  nomList.add(*funcCust) ;
1208  } else {
1209  nomList.add(*func) ;
1210  }
1211 
1212 
1213  } else {
1214 
1215  // Find the ratio term
1216  RooAbsReal* func = num;
1217  // If top level object is integral, navigate to integrand
1218  if (func->InheritsFrom(RooRealIntegral::Class())) {
1219  func = (RooAbsReal*) &((RooRealIntegral*)(func))->integrand();
1220  }
1221  if (func->InheritsFrom(RooProduct::Class())) {
1222 // cout << "product term found: " ; func->Print() ;
1223  RooArgSet comps(((RooProduct*)(func))->components()) ;
1224  RooFIter iter = comps.fwdIterator() ;
1225  RooAbsArg* arg ;
1226  while((arg=(RooAbsArg*)iter.next())) {
1227  if (arg->getAttribute("RATIO_TERM")) {
1228  ratio = (RooFormulaVar*)(arg) ;
1229  } else {
1230  origNumTerm.add(*arg) ;
1231  }
1232  }
1233  }
1234 
1235  if (ratio) {
1236 // cout << "Found ratio term in numerator: " << ratio->GetName() << endl ;
1237 // cout << "Adding only original term to numerator: " << origNumTerm << endl ;
1238  nomList.add(origNumTerm) ;
1239  } else {
1240  nomList.add(*num) ;
1241  }
1242 
1243  }
1244 
1245  for (list<string>::iterator iter = rangeComps.begin() ; iter != rangeComps.end() ; ++iter) {
1246  // If denominator is an integral, make a clone with the integration range adjusted to
1247  // the selected component of the normalization integral
1248 // cout << "NOW PROCESSING DENOMINATOR " << den->IsA()->GetName() << "::" << den->GetName() << endl ;
1249 
1250  if (string("SPECINT")==part->getStringAttribute("PROD_TERM_TYPE")) {
1251 
1252 // cout << "create integral: SPECINT case" << endl ;
1253  RooRealIntegral* orig = (RooRealIntegral*) num;
1254  RooFormulaVar* specRatio = (RooFormulaVar*) &orig->integrand() ;
1255  specIntDeps.add(orig->intVars()) ;
1256  if (orig->intRange()) {
1257  specIntRange = orig->intRange() ;
1258  }
1259  //RooProduct* numtmp = (RooProduct*) specRatio->getParameter(0) ;
1260  RooProduct* dentmp = (RooProduct*) specRatio->getParameter(1) ;
1261 
1262 // cout << "numtmp = " << numtmp->IsA()->GetName() << "::" << numtmp->GetName() << endl ;
1263 // cout << "dentmp = " << dentmp->IsA()->GetName() << "::" << dentmp->GetName() << endl ;
1264 
1265 // cout << "denominator components are " << dentmp->components() << endl ;
1266  RooArgSet comps(dentmp->components()) ;
1267  RooFIter piter = comps.fwdIterator() ;
1268  RooAbsReal* parg ;
1269  while((parg=(RooAbsReal*)piter.next())) {
1270 // cout << "now processing denominator component " << parg->IsA()->GetName() << "::" << parg->GetName() << endl ;
1271 
1272  if (ratio && parg->dependsOn(*ratio)) {
1273 // cout << "depends in value of ratio" << endl ;
1274 
1275  // Make specialize ratio instance
1276  RooAbsReal* specializedRatio = specializeRatio(*(RooFormulaVar*)ratio,iter->c_str()) ;
1277 
1278 // cout << "specRatio = " << endl ;
1279 // specializedRatio->printComponentTree() ;
1280 
1281  // Replace generic ratio with specialized ratio
1282  RooAbsArg *partCust(0) ;
1283  if (parg->InheritsFrom(RooAddition::Class())) {
1284 
1285 
1286 
1287  RooAddition* tmpadd = (RooAddition*)(parg) ;
1288 
1289  RooCustomizer cust(*tmpadd->list1().first(),Form("blah_%s",iter->c_str())) ;
1290  cust.replaceArg(*ratio,*specializedRatio) ;
1291  partCust = cust.build() ;
1292 
1293  } else {
1294  RooCustomizer cust(*parg,Form("blah_%s",iter->c_str())) ;
1295  cust.replaceArg(*ratio,*specializedRatio) ;
1296  partCust = cust.build() ;
1297  }
1298 
1299  // Print customized denominator
1300 // cout << "customized function = " << endl ;
1301 // partCust->printComponentTree() ;
1302 
1303  RooAbsReal* specializedPartCust = specializeIntegral(*(RooAbsReal*)partCust,iter->c_str()) ;
1304 
1305  // Finally divide again by ratio
1306  string name = Form("%s_divided_by_ratio",specializedPartCust->GetName()) ;
1307  RooFormulaVar* specIntFinal = new RooFormulaVar(name.c_str(),"@0/@1",RooArgList(*specializedPartCust,*specializedRatio)) ;
1308 
1309  denListList[*iter].add(*specIntFinal) ;
1310  } else {
1311 
1312 // cout << "does NOT depend on value of ratio" << endl ;
1313 // parg->Print("t") ;
1314 
1315  denListList[*iter].add(*specializeIntegral(*parg,iter->c_str())) ;
1316 
1317  }
1318  }
1319 // cout << "end iteration over denominator components" << endl ;
1320  } else {
1321 
1322  if (ratio) {
1323 
1324  RooAbsReal* specRatio = specializeRatio(*(RooFormulaVar*)ratio,iter->c_str()) ;
1325 
1326  // If integral is 'Int r(y)*g(y) dy ' then divide a posteriori by r(y)
1327 // cout << "have ratio, orig den = " << den->GetName() << endl ;
1328 
1329  RooArgSet tmp(origNumTerm) ;
1330  tmp.add(*specRatio) ;
1331  const string pname = makeRGPPName("PROD",tmp,RooArgSet(),RooArgSet(),0) ;
1332  RooProduct* specDenProd = new RooProduct(pname.c_str(),pname.c_str(),tmp) ;
1333  RooAbsReal* specInt(0) ;
1334 
1335  if (den->InheritsFrom(RooRealIntegral::Class())) {
1336  specInt = specDenProd->createIntegral(((RooRealIntegral*)den)->intVars(),iter->c_str()) ;
1337  } else if (den->InheritsFrom(RooAddition::Class())) {
1338  RooAddition* orig = (RooAddition*)den ;
1339  RooRealIntegral* origInt = (RooRealIntegral*) orig->list1().first() ;
1340  specInt = specDenProd->createIntegral(origInt->intVars(),iter->c_str()) ;
1341  } else {
1342  throw string("this should not happen") ;
1343  }
1344 
1345  //RooAbsReal* specInt = specializeIntegral(*den,iter->c_str()) ;
1346  string name = Form("%s_divided_by_ratio",specInt->GetName()) ;
1347  RooFormulaVar* specIntFinal = new RooFormulaVar(name.c_str(),"@0/@1",RooArgList(*specInt,*specRatio)) ;
1348  denListList[*iter].add(*specIntFinal) ;
1349  } else {
1350  denListList[*iter].add(*specializeIntegral(*den,iter->c_str())) ;
1351  }
1352 
1353  }
1354  }
1355 
1356  }
1357 
1358  // Do not rearrage terms if numerator and denominator are effectively empty
1359  if (nomList.getSize()==0) {
1360  return ;
1361  }
1362 
1363  string name = Form("%s_numerator",GetName()) ;
1364  // WVE FIX THIS (2)
1365 
1366  RooAbsReal* numerator = new RooProduct(name.c_str(),name.c_str(),nomList) ;
1367 
1368  RooArgSet products ;
1369 // cout << "nomList = " << nomList << endl ;
1370  for (map<string,RooArgSet>::iterator iter = denListList.begin() ; iter != denListList.end() ; ++iter) {
1371 // cout << "denList[" << iter->first << "] = " << iter->second << endl ;
1372  name = Form("%s_denominator_comp_%s",GetName(),iter->first.c_str()) ;
1373  // WVE FIX THIS (2)
1374  RooProduct* prod_comp = new RooProduct(name.c_str(),name.c_str(),iter->second) ;
1375  products.add(*prod_comp) ;
1376  }
1377  name = Form("%s_denominator_sum",GetName()) ;
1378  RooAbsReal* norm = new RooAddition(name.c_str(),name.c_str(),products) ;
1379  norm->addOwnedComponents(products) ;
1380 
1381  if (specIntDeps.getSize()>0) {
1382  // Apply posterior integration required for SPECINT case
1383 
1384  string namesr = Form("SPEC_RATIO(%s,%s)",numerator->GetName(),norm->GetName()) ;
1385  RooFormulaVar* ndr = new RooFormulaVar(namesr.c_str(),"@0/@1",RooArgList(*numerator,*norm)) ;
1386 
1387  // Integral of ratio
1388  RooAbsReal* numtmp = ndr->createIntegral(specIntDeps,specIntRange.c_str()) ;
1389 
1390  numerator = numtmp ;
1391  norm = (RooAbsReal*) RooFit::RooConst(1).Clone() ;
1392  }
1393 
1394 
1395 // cout << "numerator" << endl ;
1396 // numerator->printComponentTree("",0,5) ;
1397 // cout << "denominator" << endl ;
1398 // norm->printComponentTree("",0,5) ;
1399 
1400 
1401  // WVE DEBUG
1402  //RooMsgService::instance().debugWorkspace()->import(RooArgSet(*numerator,*norm)) ;
1403 
1404  cache._rearrangedNum.reset(numerator);
1405  cache._rearrangedDen.reset(norm);
1406  cache._isRearranged = kTRUE ;
1407 
1408 }
1409 
1410 
1411 ////////////////////////////////////////////////////////////////////////////////
1412 
1413 RooAbsReal* RooProdPdf::specializeRatio(RooFormulaVar& input, const char* targetRangeName) const
1414 {
1415  RooRealIntegral* numint = (RooRealIntegral*) input.getParameter(0) ;
1416  RooRealIntegral* denint = (RooRealIntegral*) input.getParameter(1) ;
1417 
1418  RooAbsReal* numint_spec = specializeIntegral(*numint,targetRangeName) ;
1419 
1420  RooAbsReal* ret = new RooFormulaVar(Form("ratio(%s,%s)",numint_spec->GetName(),denint->GetName()),"@0/@1",RooArgList(*numint_spec,*denint)) ;
1421  ret->addOwnedComponents(*numint_spec) ;
1422 
1423  return ret ;
1424 }
1425 
1426 
1427 
1428 ////////////////////////////////////////////////////////////////////////////////
1429 
1430 RooAbsReal* RooProdPdf::specializeIntegral(RooAbsReal& input, const char* targetRangeName) const
1431 {
1432  if (input.InheritsFrom(RooRealIntegral::Class())) {
1433 
1434  // If input is integral, recreate integral but override integration range to be targetRangeName
1435  RooRealIntegral* orig = (RooRealIntegral*)&input ;
1436 // cout << "creating integral: integrand = " << orig->integrand().GetName() << " vars = " << orig->intVars() << " range = " << targetRangeName << endl ;
1437  return orig->integrand().createIntegral(orig->intVars(),targetRangeName) ;
1438 
1439  } else if (input.InheritsFrom(RooAddition::Class())) {
1440 
1441  // If input is sum of integrals, recreate integral from first component of set, but override integration range to be targetRangeName
1442  RooAddition* orig = (RooAddition*)&input ;
1443  RooRealIntegral* origInt = (RooRealIntegral*) orig->list1().first() ;
1444 // cout << "creating integral from addition: integrand = " << origInt->integrand().GetName() << " vars = " << origInt->intVars() << " range = " << targetRangeName << endl ;
1445  return origInt->integrand().createIntegral(origInt->intVars(),targetRangeName) ;
1446 
1447  } else {
1448 
1449 // cout << "specializeIntegral: unknown input type " << input.IsA()->GetName() << "::" << input.GetName() << endl ;
1450  }
1451 
1452  return &input ;
1453 }
1454 
1455 
1456 ////////////////////////////////////////////////////////////////////////////////
1457 /// Group product into terms that can be calculated independently
1458 
1459 void RooProdPdf::groupProductTerms(std::list<std::vector<RooArgSet*>>& groupedTerms, RooArgSet& outerIntDeps,
1460  const RooLinkedList& terms, const RooLinkedList& norms,
1461  const RooLinkedList& imps, const RooLinkedList& ints, const RooLinkedList& /*cross*/) const
1462 {
1463  // Start out with each term in its own group
1464  RooFIter tIter = terms.fwdIterator() ;
1465  RooArgSet* term ;
1466  while((term=(RooArgSet*)tIter.next())) {
1467  groupedTerms.emplace_back();
1468  groupedTerms.back().emplace_back(term) ;
1469  }
1470 
1471  // Make list of imported dependents that occur in any term
1472  RooArgSet allImpDeps ;
1473  RooFIter iIter = imps.fwdIterator() ;
1474  RooArgSet *impDeps ;
1475  while((impDeps=(RooArgSet*)iIter.next())) {
1476  allImpDeps.add(*impDeps,kFALSE) ;
1477  }
1478 
1479  // Make list of integrated dependents that occur in any term
1480  RooArgSet allIntDeps ;
1481  iIter = ints.fwdIterator() ;
1482  RooArgSet *intDeps ;
1483  while((intDeps=(RooArgSet*)iIter.next())) {
1484  allIntDeps.add(*intDeps,kFALSE) ;
1485  }
1486 
1487  RooArgSet* tmp = (RooArgSet*) allIntDeps.selectCommon(allImpDeps) ;
1488  outerIntDeps.removeAll() ;
1489  outerIntDeps.add(*tmp) ;
1490  delete tmp ;
1491 
1492  // Now iteratively merge groups that should be (partially) integrated together
1493  RooFIter oidIter = outerIntDeps.fwdIterator() ;
1494  RooAbsArg* outerIntDep ;
1495  while ((outerIntDep =(RooAbsArg*)oidIter.next())) {
1496 
1497  // Collect groups that feature this dependent
1498  std::vector<RooArgSet*>* newGroup = nullptr ;
1499 
1500  // Loop over groups
1501  Bool_t needMerge = kFALSE ;
1502  auto group = groupedTerms.begin();
1503  auto nGroups = groupedTerms.size();
1504  for (size_t iGroup = 0; iGroup < nGroups; ++iGroup) {
1505 
1506  // See if any term in this group depends in any ay on outerDepInt
1507  for (auto const& term2 : *group) {
1508 
1509  Int_t termIdx = terms.IndexOf(term2) ;
1510  RooArgSet* termNormDeps = (RooArgSet*) norms.At(termIdx) ;
1511  RooArgSet* termIntDeps = (RooArgSet*) ints.At(termIdx) ;
1512  RooArgSet* termImpDeps = (RooArgSet*) imps.At(termIdx) ;
1513 
1514  if (termNormDeps->contains(*outerIntDep) ||
1515  termIntDeps->contains(*outerIntDep) ||
1516  termImpDeps->contains(*outerIntDep)) {
1517  needMerge = kTRUE ;
1518  }
1519 
1520  }
1521 
1522  if (needMerge) {
1523  // Create composite group if not yet existing
1524  if (newGroup==nullptr) {
1525  groupedTerms.emplace_back() ;
1526  newGroup = &groupedTerms.back() ;
1527  }
1528 
1529  // Add terms of this group to new term
1530  for (auto& term2 : *group) {
1531  newGroup->emplace_back(term2) ;
1532  }
1533 
1534  // Remove this group from list and delete it (but not its contents)
1535  group = groupedTerms.erase(group);
1536  } else {
1537  ++group;
1538  }
1539  }
1540 
1541  }
1542 }
1543 
1544 
1545 
1546 ////////////////////////////////////////////////////////////////////////////////
1547 /// Calculate integrals of factorized product terms over observables iset while normalized
1548 /// to observables in nset.
1549 
1550 std::vector<RooAbsReal*> RooProdPdf::processProductTerm(const RooArgSet* nset, const RooArgSet* iset, const char* isetRangeName,
1551  const RooArgSet* term,const RooArgSet& termNSet, const RooArgSet& termISet,
1552  Bool_t& isOwned, Bool_t forceWrap) const
1553 {
1554 // cout << " FOLKERT::RooProdPdf(" << GetName() <<") processProductTerm nset = " << (nset?*nset:RooArgSet()) << endl
1555 // << " _normRange = " << _normRange << endl
1556 // << " iset = " << (iset?*iset:RooArgSet()) << endl
1557 // << " isetRangeName = " << (isetRangeName?isetRangeName:"<null>") << endl
1558 // << " term = " << (term?*term:RooArgSet()) << endl
1559 // << " termNSet = " << termNSet << endl
1560 // << " termISet = " << termISet << endl
1561 // << " isOwned = " << isOwned << endl
1562 // << " forceWrap = " << forceWrap << endl ;
1563 
1564  vector<RooAbsReal*> ret(3) ; ret[0] = 0 ; ret[1] = 0 ; ret[2] = 0 ;
1565 
1566  // CASE I: factorizing term: term is integrated over all normalizing observables
1567  // -----------------------------------------------------------------------------
1568  // Check if all observbales of this term are integrated. If so the term cancels
1569  if (termNSet.getSize()>0 && termNSet.getSize()==termISet.getSize() && isetRangeName==0) {
1570 
1571 
1572  //cout << "processProductTerm(" << GetName() << ") case I " << endl ;
1573 
1574  // Term factorizes
1575  return ret ;
1576  }
1577 
1578  // CASE II: Dropped terms: if term is entirely unnormalized, it should be dropped
1579  // ------------------------------------------------------------------------------
1580  if (nset && termNSet.getSize()==0) {
1581 
1582  //cout << "processProductTerm(" << GetName() << ") case II " << endl ;
1583 
1584  // Drop terms that are not asked to be normalized
1585  return ret ;
1586  }
1587 
1588  if (iset && termISet.getSize()>0) {
1589  if (term->getSize()==1) {
1590 
1591  // CASE IIIa: Normalized and partially integrated single PDF term
1592  //---------------------------------------------------------------
1593 
1594  RooAbsPdf* pdf = (RooAbsPdf*) term->first() ;
1595 
1596  RooAbsReal* partInt = pdf->createIntegral(termISet,termNSet,isetRangeName) ;
1597  partInt->setOperMode(operMode()) ;
1598  partInt->setStringAttribute("PROD_TERM_TYPE","IIIa") ;
1599 
1600  isOwned=kTRUE ;
1601 
1602  //cout << "processProductTerm(" << GetName() << ") case IIIa func = " << partInt->GetName() << endl ;
1603 
1604  ret[0] = partInt ;
1605 
1606  // Split mode results
1607  ret[1] = pdf->createIntegral(termISet,isetRangeName) ;
1608  ret[2] = pdf->createIntegral(termNSet,normRange()) ;
1609 
1610  return ret ;
1611 
1612 
1613  } else {
1614 
1615  // CASE IIIb: Normalized and partially integrated composite PDF term
1616  //---------------------------------------------------------------
1617 
1618  // Use auxiliary class RooGenProdProj to calculate this term
1619  const std::string name = makeRGPPName("GENPROJ_",*term,termISet,termNSet,isetRangeName) ;
1620  RooAbsReal* partInt = new RooGenProdProj(name.c_str(),name.c_str(),*term,termISet,termNSet,isetRangeName) ;
1621  partInt->setStringAttribute("PROD_TERM_TYPE","IIIb") ;
1622  partInt->setOperMode(operMode()) ;
1623 
1624  //cout << "processProductTerm(" << GetName() << ") case IIIb func = " << partInt->GetName() << endl ;
1625 
1626  isOwned=kTRUE ;
1627  ret[0] = partInt ;
1628 
1629  const std::string name1 = makeRGPPName("PROD",*term,RooArgSet(),RooArgSet(),0) ;
1630 
1631  // WVE FIX THIS
1632  RooProduct* tmp_prod = new RooProduct(name1.c_str(),name1.c_str(),*term) ;
1633 
1634  ret[1] = tmp_prod->createIntegral(termISet,isetRangeName) ;
1635  ret[2] = tmp_prod->createIntegral(termNSet,normRange()) ;
1636 
1637  return ret ;
1638  }
1639  }
1640 
1641  // CASE IVa: Normalized non-integrated composite PDF term
1642  // -------------------------------------------------------
1643  if (nset && nset->getSize()>0 && term->getSize()>1) {
1644  // Composite term needs normalized integration
1645 
1646  const std::string name = makeRGPPName("GENPROJ_",*term,termISet,termNSet,isetRangeName) ;
1647  RooAbsReal* partInt = new RooGenProdProj(name.c_str(),name.c_str(),*term,termISet,termNSet,isetRangeName,normRange()) ;
1649 
1650  partInt->setStringAttribute("PROD_TERM_TYPE","IVa") ;
1651  partInt->setOperMode(operMode()) ;
1652 
1653  //cout << "processProductTerm(" << GetName() << ") case IVa func = " << partInt->GetName() << endl ;
1654 
1655  isOwned=kTRUE ;
1656  ret[0] = partInt ;
1657 
1658  const std::string name1 = makeRGPPName("PROD",*term,RooArgSet(),RooArgSet(),0) ;
1659 
1660  // WVE FIX THIS
1661  RooProduct* tmp_prod = new RooProduct(name1.c_str(),name1.c_str(),*term) ;
1662 
1663  ret[1] = tmp_prod->createIntegral(termISet,isetRangeName) ;
1664  ret[2] = tmp_prod->createIntegral(termNSet,normRange()) ;
1665 
1666  return ret ;
1667  }
1668 
1669  // CASE IVb: Normalized, non-integrated single PDF term
1670  // -----------------------------------------------------
1671  RooFIter pIter = term->fwdIterator() ;
1672  RooAbsPdf* pdf ;
1673  while((pdf=(RooAbsPdf*)pIter.next())) {
1674 
1675  if (forceWrap) {
1676 
1677  // Construct representative name of normalization wrapper
1678  TString name(pdf->GetName()) ;
1679  name.Append("_NORM[") ;
1680  RooFIter nIter = termNSet.fwdIterator() ;
1681  RooAbsArg* arg ;
1682  Bool_t first(kTRUE) ;
1683  while((arg=(RooAbsArg*)nIter.next())) {
1684  if (!first) {
1685  name.Append(",") ;
1686  } else {
1687  first=kFALSE ;
1688  }
1689  name.Append(arg->GetName()) ;
1690  }
1691  if (normRange()) {
1692  name.Append("|") ;
1693  name.Append(normRange()) ;
1694  }
1695  name.Append("]") ;
1696 
1697  RooAbsReal* partInt = new RooRealIntegral(name.Data(),name.Data(),*pdf,RooArgSet(),&termNSet) ;
1698  partInt->setStringAttribute("PROD_TERM_TYPE","IVb") ;
1699  isOwned=kTRUE ;
1700 
1701  //cout << "processProductTerm(" << GetName() << ") case IVb func = " << partInt->GetName() << endl ;
1702 
1703  ret[0] = partInt ;
1704 
1705  ret[1] = pdf->createIntegral(RooArgSet()) ;
1706  ret[2] = pdf->createIntegral(termNSet,normRange()) ;
1707 
1708  return ret ;
1709 
1710 
1711  } else {
1712  isOwned=kFALSE ;
1713 
1714  //cout << "processProductTerm(" << GetName() << ") case IVb func = " << pdf->GetName() << endl ;
1715 
1716 
1717  pdf->setStringAttribute("PROD_TERM_TYPE","IVb") ;
1718  ret[0] = pdf ;
1719 
1720  ret[1] = pdf->createIntegral(RooArgSet()) ;
1721  ret[2] = termNSet.getSize()>0 ? pdf->createIntegral(termNSet,normRange()) : ((RooAbsReal*)RooFit::RooConst(1).clone("1")) ;
1722  return ret ;
1723  }
1724  }
1725 
1726  coutE(Eval) << "RooProdPdf::processProductTerm(" << GetName() << ") unidentified term!!!" << endl ;
1727  return ret ;
1728 }
1729 
1730 
1731 
1732 
1733 ////////////////////////////////////////////////////////////////////////////////
1734 /// Make an appropriate automatic name for a RooGenProdProj object in getPartIntList()
1735 
1736 std::string RooProdPdf::makeRGPPName(const char* pfx, const RooArgSet& term, const RooArgSet& iset,
1737  const RooArgSet& nset, const char* isetRangeName) const
1738 {
1739  // Make an appropriate automatic name for a RooGenProdProj object in getPartIntList()
1740 
1741  std::ostringstream os(pfx);
1742  os << "[";
1743 
1744  RooFIter pIter = term.fwdIterator() ;
1745  // Encode component names
1746  Bool_t first(kTRUE) ;
1747  RooAbsPdf* pdf ;
1748  while ((pdf=(RooAbsPdf*)pIter.next())) {
1749  if (!first) os << "_X_";
1750  first = kFALSE;
1751  os << pdf->GetName();
1752  }
1753  os << "]" << integralNameSuffix(iset,&nset,isetRangeName,kTRUE);
1754 
1755  return os.str();
1756 }
1757 
1758 
1759 
1760 ////////////////////////////////////////////////////////////////////////////////
1761 /// Force RooRealIntegral to offer all observables for internal integration
1762 
1764 {
1765  return kTRUE ;
1766 }
1767 
1768 
1769 
1770 ////////////////////////////////////////////////////////////////////////////////
1771 /// Determine which part (if any) of given integral can be performed analytically.
1772 /// If any analytical integration is possible, return integration scenario code.
1773 ///
1774 /// RooProdPdf implements two strategies in implementing analytical integrals
1775 ///
1776 /// First, PDF components whose entire set of dependents are requested to be integrated
1777 /// can be dropped from the product, as they will integrate out to 1 by construction
1778 ///
1779 /// Second, RooProdPdf queries each remaining component PDF for its analytical integration
1780 /// capability of the requested set ('allVars'). It finds the largest common set of variables
1781 /// that can be integrated by all remaining components. If such a set exists, it reconfirms that
1782 /// each component is capable of analytically integrating the common set, and combines the components
1783 /// individual integration codes into a single integration code valid for RooProdPdf.
1784 
1786  const RooArgSet* normSet, const char* rangeName) const
1787 {
1788  if (_forceNumInt) return 0 ;
1789 
1790  // Declare that we can analytically integrate all requested observables
1791  analVars.add(allVars) ;
1792 
1793  // Retrieve (or create) the required partial integral list
1794  Int_t code = getPartIntList(normSet,&allVars,rangeName);
1795 
1796  return code+1 ;
1797 }
1798 
1799 
1800 
1801 
1802 ////////////////////////////////////////////////////////////////////////////////
1803 /// Return analytical integral defined by given scenario code
1804 
1805 Double_t RooProdPdf::analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* rangeName) const
1806 {
1807  // No integration scenario
1808  if (code==0) {
1809  return getVal(normSet) ;
1810  }
1811 
1812 
1813  // WVE needs adaptation for rangename feature
1814 
1815  // Partial integration scenarios
1816  CacheElem* cache = (CacheElem*) _cacheMgr.getObjByIndex(code-1) ;
1817 
1818  // If cache has been sterilized, revive this slot
1819  if (cache==0) {
1820  std::unique_ptr<RooArgSet> vars{getParameters(RooArgSet())} ;
1821  RooArgSet nset = _cacheMgr.selectFromSet1(*vars, code-1) ;
1822  RooArgSet iset = _cacheMgr.selectFromSet2(*vars, code-1) ;
1823 
1824  Int_t code2 = getPartIntList(&nset, &iset, rangeName) ;
1825 
1826  // preceding call to getPartIntList guarantees non-null return
1827  // coverity[NULL_RETURNS]
1828  cache = (CacheElem*) _cacheMgr.getObj(&nset,&iset,&code2,rangeName) ;
1829  }
1830 
1831  Double_t val = calculate(*cache,kTRUE) ;
1832 // cout << "RPP::aIWN(" << GetName() << ") ,code = " << code << ", value = " << val << endl ;
1833 
1834  return val ;
1835 }
1836 
1837 
1838 
1839 ////////////////////////////////////////////////////////////////////////////////
1840 /// Obsolete
1841 
1843 { return kFALSE ; }
1844 
1845 
1846 
1847 
1848 ////////////////////////////////////////////////////////////////////////////////
1849 /// If this product contains exactly one extendable p.d.f return the extension abilities of
1850 /// that p.d.f, otherwise return CanNotBeExtended
1851 
1853 {
1855 }
1856 
1857 
1858 
1859 ////////////////////////////////////////////////////////////////////////////////
1860 /// Return the expected number of events associated with the extendable input PDF
1861 /// in the product. If there is no extendable term, abort.
1862 
1864 {
1865  if (_extendedIndex<0) {
1866  coutF(Generation) << "Requesting expected number of events from a RooProdPdf that does not contain an extended p.d.f" << endl ;
1867  throw std::logic_error(std::string("RooProdPdf ") + GetName() + " could not be extended.");
1868  }
1869 
1870  return ((RooAbsPdf*)_pdfList.at(_extendedIndex))->expectedEvents(nset) ;
1871 }
1872 
1873 
1874 
1875 ////////////////////////////////////////////////////////////////////////////////
1876 /// Return generator context optimized for generating events from product p.d.f.s
1877 
1879  const RooArgSet* auxProto, Bool_t verbose) const
1880 {
1881  if (_useDefaultGen) return RooAbsPdf::genContext(vars,prototype,auxProto,verbose) ;
1882  return new RooProdGenContext(*this,vars,prototype,auxProto,verbose) ;
1883 }
1884 
1885 
1886 
1887 ////////////////////////////////////////////////////////////////////////////////
1888 /// Query internal generation capabilities of component p.d.f.s and aggregate capabilities
1889 /// into master configuration passed to the generator context
1890 
1891 Int_t RooProdPdf::getGenerator(const RooArgSet& directVars, RooArgSet &generateVars, Bool_t staticInitOK) const
1892 {
1893  if (!_useDefaultGen) return 0 ;
1894 
1895  // Find the subset directVars that only depend on a single PDF in the product
1896  RooArgSet directSafe ;
1897  RooFIter dIter = directVars.fwdIterator() ;
1898  RooAbsArg* arg ;
1899  while((arg=(RooAbsArg*)dIter.next())) {
1900  if (isDirectGenSafe(*arg)) directSafe.add(*arg) ;
1901  }
1902 
1903 
1904  // Now find direct integrator for relevant components ;
1905  RooAbsPdf* pdf ;
1906  std::vector<Int_t> code;
1907  code.reserve(64);
1908  RooFIter pdfIter = _pdfList.fwdIterator();
1909  while((pdf=(RooAbsPdf*)pdfIter.next())) {
1910  RooArgSet pdfDirect ;
1911  Int_t pdfCode = pdf->getGenerator(directSafe,pdfDirect,staticInitOK);
1912  code.push_back(pdfCode);
1913  if (pdfCode != 0) {
1914  generateVars.add(pdfDirect) ;
1915  }
1916  }
1917 
1918 
1919  if (generateVars.getSize()>0) {
1920  Int_t masterCode = _genCode.store(code) ;
1921  return masterCode+1 ;
1922  } else {
1923  return 0 ;
1924  }
1925 }
1926 
1927 
1928 
1929 ////////////////////////////////////////////////////////////////////////////////
1930 /// Forward one-time initialization call to component generation initialization
1931 /// methods.
1932 
1934 {
1935  if (!_useDefaultGen) return ;
1936 
1937  const std::vector<Int_t>& codeList = _genCode.retrieve(code-1) ;
1938  RooAbsPdf* pdf ;
1939  Int_t i(0) ;
1940  RooFIter pdfIter = _pdfList.fwdIterator();
1941  while((pdf=(RooAbsPdf*)pdfIter.next())) {
1942  if (codeList[i]!=0) {
1943  pdf->initGenerator(codeList[i]) ;
1944  }
1945  i++ ;
1946  }
1947 }
1948 
1949 
1950 
1951 ////////////////////////////////////////////////////////////////////////////////
1952 /// Generate a single event with configuration specified by 'code'
1953 /// Defer internal generation to components as encoded in the _genCode
1954 /// registry for given generator code.
1955 
1957 {
1958  if (!_useDefaultGen) return ;
1959 
1960  const std::vector<Int_t>& codeList = _genCode.retrieve(code-1) ;
1961  RooAbsPdf* pdf ;
1962  Int_t i(0) ;
1963  RooFIter pdfIter = _pdfList.fwdIterator();
1964  while((pdf=(RooAbsPdf*)pdfIter.next())) {
1965  if (codeList[i]!=0) {
1966  pdf->generateEvent(codeList[i]) ;
1967  }
1968  i++ ;
1969  }
1970 }
1971 
1972 
1973 
1974 ////////////////////////////////////////////////////////////////////////////////
1975 /// Return RooAbsArg components contained in the cache
1976 
1978 {
1979  RooArgList ret ;
1980  ret.add(_partList) ;
1981  ret.add(_numList) ;
1982  ret.add(_denList) ;
1983  if (_rearrangedNum) ret.add(*_rearrangedNum) ;
1984  if (_rearrangedDen) ret.add(*_rearrangedDen) ;
1985  return ret ;
1986 
1987 }
1988 
1989 
1990 
1991 ////////////////////////////////////////////////////////////////////////////////
1992 /// Hook function to print cache contents in tree printing of RooProdPdf
1993 
1994 void RooProdPdf::CacheElem::printCompactTreeHook(ostream& os, const char* indent, Int_t curElem, Int_t maxElem)
1995 {
1996  if (curElem==0) {
1997  os << indent << "RooProdPdf begin partial integral cache" << endl ;
1998  }
1999 
2000  auto indent2 = std::string(indent) + "[" + std::to_string(curElem) + "]";
2001  for(auto const& arg : _partList) {
2002  arg->printCompactTree(os,indent2.c_str()) ;
2003  }
2004 
2005  if (curElem==maxElem) {
2006  os << indent << "RooProdPdf end partial integral cache" << endl ;
2007  }
2008 }
2009 
2010 
2011 
2012 ////////////////////////////////////////////////////////////////////////////////
2013 /// Forward determination of safety of internal generator code to
2014 /// component p.d.f that would generate the given observable
2015 
2017 {
2018  // Only override base class behaviour if default generator method is enabled
2019  if (!_useDefaultGen) return RooAbsPdf::isDirectGenSafe(arg) ;
2020 
2021  // Argument may appear in only one PDF component
2022  RooAbsPdf* pdf, *thePdf(0) ;
2023  RooFIter pdfIter = _pdfList.fwdIterator();
2024  while((pdf=(RooAbsPdf*)pdfIter.next())) {
2025 
2026  if (pdf->dependsOn(arg)) {
2027  // Found PDF depending on arg
2028 
2029  // If multiple PDFs depend on arg directGen is not safe
2030  if (thePdf) return kFALSE ;
2031 
2032  thePdf = pdf ;
2033  }
2034  }
2035  // Forward call to relevant component PDF
2036  return thePdf?(thePdf->isDirectGenSafe(arg)):kFALSE ;
2037 }
2038 
2039 
2040 
2041 ////////////////////////////////////////////////////////////////////////////////
2042 /// Look up user specified normalization set for given input PDF component
2043 
2045 {
2046  Int_t idx = _pdfList.index(&pdf) ;
2047  if (idx<0) return 0 ;
2048  return (RooArgSet*) _pdfNSetList.At(idx) ;
2049 }
2050 
2051 
2052 
2053 ////////////////////////////////////////////////////////////////////////////////
2054 /// Return all parameter constraint p.d.f.s on parameters listed in constrainedParams.
2055 /// The observables set is required to distinguish unambiguously p.d.f in terms
2056 /// of observables and parameters, which are not constraints, and p.d.fs in terms
2057 /// of parameters only, which can serve as constraints p.d.f.s
2058 
2059 RooArgSet* RooProdPdf::getConstraints(const RooArgSet& observables, RooArgSet& constrainedParams, Bool_t stripDisconnected) const
2060 {
2061  RooArgSet constraints ;
2062  RooArgSet pdfParams, conParams ;
2063 
2064  // Loop over p.d.f. components
2065  RooFIter piter = _pdfList.fwdIterator() ;
2066  RooAbsPdf* pdf ;
2067  while((pdf=(RooAbsPdf*)piter.next())) {
2068  // A constraint term is a p.d.f that does not depend on any of the listed observables
2069  // but does depends on any of the parameters that should be constrained
2070  if (!pdf->dependsOnValue(observables) && pdf->dependsOnValue(constrainedParams)) {
2071  constraints.add(*pdf) ;
2072  RooArgSet* tmp = pdf->getParameters(observables) ;
2073  conParams.add(*tmp,kTRUE) ;
2074  delete tmp ;
2075  } else {
2076  RooArgSet* tmp = pdf->getParameters(observables) ;
2077  pdfParams.add(*tmp,kTRUE) ;
2078  delete tmp ;
2079  }
2080  }
2081 
2082  // Strip any constraints that are completely decoupled from the other product terms
2083  RooArgSet* finalConstraints = new RooArgSet("constraints") ;
2084  RooFIter citer = constraints.fwdIterator() ;
2085  while((pdf=(RooAbsPdf*)citer.next())) {
2086  if (pdf->dependsOnValue(pdfParams) || !stripDisconnected) {
2087  finalConstraints->add(*pdf) ;
2088  } else {
2089  coutI(Minimization) << "RooProdPdf::getConstraints(" << GetName() << ") omitting term " << pdf->GetName()
2090  << " as constraint term as it does not share any parameters with the other pdfs in product. "
2091  << "To force inclusion in likelihood, add an explicit Constrain() argument for the target parameter" << endl ;
2092  }
2093  }
2094 
2095  // Now remove from constrainedParams all parameters that occur exclusively in constraint term and not in regular pdf term
2096 
2097  RooArgSet* cexl = (RooArgSet*) conParams.selectCommon(constrainedParams) ;
2098  cexl->remove(pdfParams,kTRUE,kTRUE) ;
2099  constrainedParams.remove(*cexl,kTRUE,kTRUE) ;
2100  delete cexl ;
2101 
2102  return finalConstraints ;
2103 }
2104 
2105 
2106 
2107 
2108 ////////////////////////////////////////////////////////////////////////////////
2109 /// Return all parameter constraint p.d.f.s on parameters listed in constrainedParams.
2110 /// The observables set is required to distinguish unambiguously p.d.f in terms
2111 /// of observables and parameters, which are not constraints, and p.d.fs in terms
2112 /// of parameters only, which can serve as constraints p.d.f.s
2113 
2115 {
2116  RooArgSet* connectedPars = new RooArgSet("connectedPars") ;
2117  for (const auto arg : _pdfList) {
2118  // Check if term is relevant
2119  if (arg->dependsOn(observables)) {
2120  RooArgSet* tmp = arg->getParameters(observables) ;
2121  connectedPars->add(*tmp) ;
2122  delete tmp ;
2123  }
2124  }
2125  return connectedPars ;
2126 }
2127 
2128 
2129 
2130 
2131 ////////////////////////////////////////////////////////////////////////////////
2132 
2133 void RooProdPdf::getParametersHook(const RooArgSet* nset, RooArgSet* params, Bool_t stripDisconnected) const
2134 {
2135  if (!stripDisconnected) return ;
2136  if (!nset || nset->getSize()==0) return ;
2137 
2138  // Get/create appropriate term list for this normalization set
2139  Int_t code = getPartIntList(nset, nullptr);
2140  RooArgList & plist = static_cast<CacheElem*>(_cacheMgr.getObj(nset, &code))->_partList;
2141 
2142  // Strip any terms from params that do not depend on any term
2143  RooArgSet tostrip ;
2144  for (auto param : *params) {
2145  Bool_t anyDep(kFALSE) ;
2146  for (auto term : plist) {
2147  if (term->dependsOnValue(*param)) {
2148  anyDep=kTRUE ;
2149  }
2150  }
2151  if (!anyDep) {
2152  tostrip.add(*param) ;
2153  }
2154  }
2155 
2156  if (tostrip.getSize()>0) {
2157  params->remove(tostrip,kTRUE,kTRUE);
2158  }
2159 
2160 }
2161 
2162 
2163 
2164 ////////////////////////////////////////////////////////////////////////////////
2165 /// Interface function used by test statistics to freeze choice of range
2166 /// for interpretation of conditional product terms
2167 
2168 void RooProdPdf::selectNormalizationRange(const char* rangeName, Bool_t force)
2169 {
2170  if (!force && _refRangeName) {
2171  return ;
2172  }
2173 
2174  fixRefRange(rangeName) ;
2175 }
2176 
2177 
2178 
2179 
2180 ////////////////////////////////////////////////////////////////////////////////
2181 
2182 void RooProdPdf::fixRefRange(const char* rangeName)
2183 {
2184  _refRangeName = (TNamed*)RooNameReg::ptr(rangeName) ;
2185 }
2186 
2187 
2188 
2189 ////////////////////////////////////////////////////////////////////////////////
2190 /// Forward the plot sampling hint from the p.d.f. that defines the observable obs
2191 
2192 std::list<Double_t>* RooProdPdf::plotSamplingHint(RooAbsRealLValue& obs, Double_t xlo, Double_t xhi) const
2193 {
2194  RooAbsPdf* pdf ;
2195  RooFIter pdfIter = _pdfList.fwdIterator();
2196  while((pdf=(RooAbsPdf*)pdfIter.next())) {
2197  list<Double_t>* hint = pdf->plotSamplingHint(obs,xlo,xhi) ;
2198  if (hint) {
2199  return hint ;
2200  }
2201  }
2202 
2203  return 0 ;
2204 }
2205 
2206 
2207 
2208 ////////////////////////////////////////////////////////////////////////////////
2209 /// If all components that depend on obs are binned that so is the product
2210 
2212 {
2213  RooAbsPdf* pdf ;
2214  RooFIter pdfIter = _pdfList.fwdIterator();
2215  while((pdf=(RooAbsPdf*)pdfIter.next())) {
2216  if (pdf->dependsOn(obs) && !pdf->isBinnedDistribution(obs)) {
2217  return kFALSE ;
2218  }
2219  }
2220 
2221  return kTRUE ;
2222 }
2223 
2224 
2225 
2226 
2227 
2228 
2229 ////////////////////////////////////////////////////////////////////////////////
2230 /// Forward the plot sampling hint from the p.d.f. that defines the observable obs
2231 
2232 std::list<Double_t>* RooProdPdf::binBoundaries(RooAbsRealLValue& obs, Double_t xlo, Double_t xhi) const
2233 {
2234  RooAbsPdf* pdf ;
2235  RooFIter pdfIter = _pdfList.fwdIterator();
2236  while((pdf=(RooAbsPdf*)pdfIter.next())) {
2237  list<Double_t>* hint = pdf->binBoundaries(obs,xlo,xhi) ;
2238  if (hint) {
2239  return hint ;
2240  }
2241  }
2242 
2243  return 0 ;
2244 }
2245 
2246 
2247 ////////////////////////////////////////////////////////////////////////////////
2248 /// Label OK'ed components of a RooProdPdf with cache-and-track, _and_ label all RooProdPdf
2249 /// descendants with extra information about (conditional) normalization, needed to be able
2250 /// to Cache-And-Track them outside the RooProdPdf context.
2251 
2253 {
2254  for (const auto parg : _pdfList) {
2255 
2256  if (parg->canNodeBeCached()==Always) {
2257  trackNodes.add(*parg) ;
2258 // cout << "tracking node RooProdPdf component " << parg << " " << parg->IsA()->GetName() << "::" << parg->GetName() << endl ;
2259 
2260  // Additional processing to fix normalization sets in case product defines conditional observables
2261  RooArgSet* pdf_nset = findPdfNSet((RooAbsPdf&)(*parg)) ;
2262  if (pdf_nset) {
2263  // Check if conditional normalization is specified
2265  if (string("nset")==pdf_nset->GetName() && pdf_nset->getSize()>0) {
2266  parg->setStringAttribute("CATNormSet",getColonSeparatedNameString(*pdf_nset).c_str()) ;
2267  }
2268  if (string("cset")==pdf_nset->GetName()) {
2269  parg->setStringAttribute("CATCondSet",getColonSeparatedNameString(*pdf_nset).c_str()) ;
2270  }
2271  } else {
2272  coutW(Optimization) << "RooProdPdf::setCacheAndTrackHints(" << GetName() << ") WARNING product pdf does not specify a normalization set for component " << parg->GetName() << endl ;
2273  }
2274  }
2275  }
2276 }
2277 
2278 
2279 
2280 ////////////////////////////////////////////////////////////////////////////////
2281 /// Customized printing of arguments of a RooProdPdf to more intuitively reflect the contents of the
2282 /// product operator construction
2283 
2284 void RooProdPdf::printMetaArgs(ostream& os) const
2285 {
2286  RooFIter niter = _pdfNSetList.fwdIterator() ;
2287  for (int i=0 ; i<_pdfList.getSize() ; i++) {
2288  if (i>0) os << " * " ;
2289  RooArgSet* ncset = (RooArgSet*) niter.next() ;
2290  os << _pdfList.at(i)->GetName() ;
2291  if (ncset->getSize()>0) {
2292  if (string("nset")==ncset->GetName()) {
2293  os << *ncset ;
2294  } else {
2295  os << "|" ;
2296  RooFIter nciter = ncset->fwdIterator() ;
2297  RooAbsArg* arg ;
2298  Bool_t first(kTRUE) ;
2299  while((arg=(RooAbsArg*)nciter.next())) {
2300  if (!first) {
2301  os << "," ;
2302  } else {
2303  first = kFALSE ;
2304  }
2305  os << arg->GetName() ;
2306  }
2307  }
2308  }
2309  }
2310  os << " " ;
2311 }
2312 
2313 
2314 
2315 ////////////////////////////////////////////////////////////////////////////////
2316 /// Implement support for node removal
2317 
2318 Bool_t RooProdPdf::redirectServersHook(const RooAbsCollection& /*newServerList*/, Bool_t /*mustReplaceAll*/, Bool_t nameChange, Bool_t /*isRecursive*/)
2319 {
2320  if (nameChange && _pdfList.find("REMOVAL_DUMMY")) {
2321 
2322  cxcoutD(LinkStateMgmt) << "RooProdPdf::redirectServersHook(" << GetName() << "): removing REMOVAL_DUMMY" << endl ;
2323 
2324  // Remove node from _pdfList proxy and remove corresponding entry from normset list
2325  RooAbsArg* pdfDel = _pdfList.find("REMOVAL_DUMMY") ;
2326 
2327  TObject* setDel = _pdfNSetList.At(_pdfList.index("REMOVAL_DUMMY")) ;
2328  _pdfList.remove(*pdfDel) ;
2329  _pdfNSetList.Remove(setDel) ;
2330 
2331  // Clear caches
2332  _cacheMgr.reset() ;
2333  }
2334  return kFALSE ;
2335 }
2336 
2337 void RooProdPdf::CacheElem::writeToStream(std::ostream& os) const {
2338  using namespace RooHelpers;
2339  os << "_partList\n";
2340  os << getColonSeparatedNameString(_partList) << "\n";
2341  os << "_numList\n";
2342  os << getColonSeparatedNameString(_numList) << "\n";
2343  os << "_denList\n";
2344  os << getColonSeparatedNameString(_denList) << "\n";
2345  os << "_ownedList\n";
2346  os << getColonSeparatedNameString(_ownedList) << "\n";
2347  os << "_normList\n";
2348  for(auto const& set : _normList) {
2349  os << getColonSeparatedNameString(*set) << "\n";
2350  }
2351  os << "_isRearranged" << "\n";
2352  os << _isRearranged << "\n";
2353  os << "_rearrangedNum" << "\n";
2354  if(_rearrangedNum) {
2355  os << getColonSeparatedNameString(*_rearrangedNum) << "\n";
2356  } else {
2357  os << "nullptr" << "\n";
2358  }
2359  os << "_rearrangedDen" << "\n";
2360  if(_rearrangedDen) {
2361  os << getColonSeparatedNameString(*_rearrangedDen) << "\n";
2362  } else {
2363  os << "nullptr" << "\n";
2364  }
2365 }
2366 
2367 void RooProdPdf::writeCacheToStream(std::ostream& os, RooArgSet const* nset) const {
2368  getCacheElem(nset)->writeToStream(os);
2369 }
RooFormulaVar.h
RooAbsArg::Clone
virtual TObject * Clone(const char *newname=0) const
Make a clone of an object using the Streamer facility.
Definition: RooAbsArg.h:84
RooCacheManager::reset
void reset()
Clear the cache.
Definition: RooCacheManager.h:193
RooProdPdf::_pdfList
RooListProxy _pdfList
Definition: RooProdPdf.h:172
l
auto * l
Definition: textangle.C:4
RooAbsPdf::_normSet
RooArgSet * _normSet
Normalization integral (owned by _normMgr)
Definition: RooAbsPdf.h:325
RooCacheManager::setObj
Int_t setObj(const RooArgSet *nset, T *obj, const TNamed *isetRangeName=0)
Definition: RooCacheManager.h:50
RooBatchCompute.h
RooAbsPdf::CacheElem
friend class CacheElem
Definition: RooAbsPdf.h:336
RooAbsPdf::_normRange
TString _normRange
MC generator configuration specific for this object.
Definition: RooAbsPdf.h:357
RooWorkspace.h
RooCmdArg
RooCmdArg is a named container for two doubles, two integers two object points and three string point...
Definition: RooCmdArg.h:27
RooAbsPdf::genContext
virtual RooAbsGenContext * genContext(const RooArgSet &vars, const RooDataSet *prototype=0, const RooArgSet *auxProto=0, Bool_t verbose=kFALSE) const
Interface function to create a generator context from a p.d.f.
Definition: RooAbsPdf.cxx:2066
RooProdGenContext.h
fit1_py.fill
fill
Definition: fit1_py.py:6
first
Definition: first.py:1
RooHelpers.h
RooFit::Minimization
@ Minimization
Definition: RooGlobalFunc.h:60
RooAbsArg::operMode
OperMode operMode() const
Query the operation mode of this node.
Definition: RooAbsArg.h:484
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:100
RooProdPdf::isBinnedDistribution
Bool_t isBinnedDistribution(const RooArgSet &obs) const
If all components that depend on obs are binned that so is the product.
Definition: RooProdPdf.cxx:2211
RooBatchCompute::RunContext::makeBatch
RooSpan< double > makeBatch(const RooAbsReal *owner, std::size_t size)
Create a writable batch.
Definition: RunContext.cxx:87
RooCmdArg::getSet
const RooArgSet * getSet(Int_t idx) const
Return RooArgSet stored in slot idx.
Definition: RooCmdArg.cxx:206
RooProdPdf::_genCode
RooAICRegistry _genCode
Definition: RooProdPdf.h:169
RooMsgService.h
RooAbsCollection::first
RooAbsArg * first() const
Definition: RooAbsCollection.h:236
RooAbsReal::plotSamplingHint
virtual std::list< Double_t > * plotSamplingHint(RooAbsRealLValue &obs, Double_t xlo, Double_t xhi) const
Interface for returning an optional hint for initial sampling points when constructing a curve projec...
Definition: RooAbsReal.cxx:3824
RooNameReg::ptr
static const TNamed * ptr(const char *stringPtr)
Return a unique TNamed pointer for given C++ string.
Definition: RooNameReg.cxx:93
RooCmdArg::getInt
Int_t getInt(Int_t idx) const
Definition: RooCmdArg.h:80
RooProdPdf::evaluate
Double_t evaluate() const
Calculate current value of object.
Definition: RooProdPdf.cxx:464
RooProdPdf::writeCacheToStream
void writeCacheToStream(std::ostream &os, RooArgSet const *nset) const
Definition: RooProdPdf.cxx:2367
RooAbsRealLValue::getMax
virtual Double_t getMax(const char *name=0) const
Get maximum of currently defined range.
Definition: RooAbsRealLValue.h:89
RooProdPdf::specializeIntegral
RooAbsReal * specializeIntegral(RooAbsReal &orig, const char *targetRangeName) const
Definition: RooProdPdf.cxx:1430
RooAbsPdf::generateEvent
virtual void generateEvent(Int_t code)
Interface for generation of an event using the algorithm corresponding to the specified code.
Definition: RooAbsPdf.cxx:2496
RooProdPdf::CacheElem::_normList
std::vector< std::unique_ptr< RooArgSet > > _normList
Definition: RooProdPdf.h:146
RooProdPdf::isDirectGenSafe
virtual Bool_t isDirectGenSafe(const RooAbsArg &arg) const
Forward determination of safety of internal generator code to component p.d.f that would generate the...
Definition: RooProdPdf.cxx:2016
RooFit::InputArguments
@ InputArguments
Definition: RooGlobalFunc.h:61
RooProdPdf::forceAnalyticalInt
virtual Bool_t forceAnalyticalInt(const RooAbsArg &dep) const
Force RooRealIntegral to offer all observables for internal integration.
Definition: RooProdPdf.cxx:1763
RooAbsArg::getComponents
RooArgSet * getComponents() const
Create a RooArgSet with all components (branch nodes) of the expression tree headed by this object.
Definition: RooAbsArg.cxx:746
RooAbsArg::Always
@ Always
Definition: RooAbsArg.h:389
RooAbsArg::getParameters
RooArgSet * getParameters(const RooAbsData *data, bool stripDisconnected=true) const
Create a list of leaf nodes in the arg tree starting with ourself as top node that don't match any of...
Definition: RooAbsArg.cxx:578
RooCacheManager::getObj
T * getObj(const RooArgSet *nset, Int_t *sterileIndex=0, const TNamed *isetRangeName=0)
Definition: RooCacheManager.h:45
RooProdPdf::_cutOff
Double_t _cutOff
Registry of composite direct generator codes.
Definition: RooProdPdf.h:171
RooAbsReal::isBinnedDistribution
virtual Bool_t isBinnedDistribution(const RooArgSet &) const
Tests if the distribution is binned. Unless overridden by derived classes, this always returns false.
Definition: RooAbsReal.h:341
TString::Data
const char * Data() const
Definition: TString.h:369
RooProdPdf::_cacheMgr
RooObjCacheManager _cacheMgr
Definition: RooProdPdf.h:155
RooAbsReal::createIntegral
RooAbsReal * createIntegral(const RooArgSet &iset, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Create an object that represents the integral of the function over one or more observables listed in ...
Definition: RooAbsReal.cxx:548
RooProdPdf::CacheElem::printCompactTreeHook
void printCompactTreeHook(std::ostream &, const char *, Int_t, Int_t) override
Hook function to print cache contents in tree printing of RooProdPdf.
Definition: RooProdPdf.cxx:1994
ClassImp
#define ClassImp(name)
Definition: Rtypes.h:364
Form
char * Form(const char *fmt,...)
RooAbsArg::dependsOnValue
Bool_t dependsOnValue(const RooAbsCollection &serverList, const RooAbsArg *ignoreArg=0) const
Check whether this object depends on values from an element in the serverList.
Definition: RooAbsArg.h:102
RooProdPdf::selectNormalizationRange
virtual void selectNormalizationRange(const char *rangeName=0, Bool_t force=kFALSE)
Interface function used by test statistics to freeze choice of range for interpretation of conditiona...
Definition: RooProdPdf.cxx:2168
RooConstVar::clone
virtual TObject * clone(const char *newname) const
Definition: RooConstVar.h:32
RooCustomizer
RooCustomizer is a factory class to produce clones of a prototype composite PDF object with the same ...
Definition: RooCustomizer.h:35
RooAbsCollection::fwdIterator
RooFIter fwdIterator() const
One-time forward iterator.
Definition: RooAbsCollection.h:193
coutE
#define coutE(a)
Definition: RooMsgService.h:33
RooAbsReal::binBoundaries
virtual std::list< Double_t > * binBoundaries(RooAbsRealLValue &obs, Double_t xlo, Double_t xhi) const
Retrieve bin boundaries if this distribution is binned in obs.
Definition: RooAbsReal.cxx:3813
coutW
#define coutW(a)
Definition: RooMsgService.h:32
RooArgList
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:21
RooProdPdf::getConstraints
virtual RooArgSet * getConstraints(const RooArgSet &observables, RooArgSet &constrainedParams, Bool_t stripDisconnected) const
Return all parameter constraint p.d.f.s on parameters listed in constrainedParams.
Definition: RooProdPdf.cxx:2059
RooAbsReal::getVal
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:91
RooFit::Optimization
@ Optimization
Definition: RooGlobalFunc.h:61
RooProdPdf::getPartIntList
Int_t getPartIntList(const RooArgSet *nset, const RooArgSet *iset, const char *isetRangeName=0) const
Return list of (partial) integrals of product terms for integration of p.d.f over observables iset wh...
Definition: RooProdPdf.cxx:782
RooProdPdf::getConnectedParameters
RooArgSet * getConnectedParameters(const RooArgSet &observables) const
Return all parameter constraint p.d.f.s on parameters listed in constrainedParams.
Definition: RooProdPdf.cxx:2114
RooProdPdf::CacheElem::_rearrangedNum
std::unique_ptr< RooAbsReal > _rearrangedNum
Definition: RooProdPdf.h:148
RooProdPdf::plotSamplingHint
virtual std::list< Double_t > * plotSamplingHint(RooAbsRealLValue &obs, Double_t xlo, Double_t xhi) const
Forward the plot sampling hint from the p.d.f. that defines the observable obs.
Definition: RooProdPdf.cxx:2192
RooAbsCollection::remove
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE)
Remove the specified argument from our list.
Definition: RooAbsCollection.cxx:585
RooAbsCollection::find
RooAbsArg * find(const char *name) const
Find object with given name in list.
Definition: RooAbsCollection.cxx:810
RooProdPdf::_refRangeName
TNamed * _refRangeName
Definition: RooProdPdf.h:179
TString::Contains
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:624
RooAbsArg::setOperMode
void setOperMode(OperMode mode, Bool_t recurseADirty=kTRUE)
Set the operation mode of this node.
Definition: RooAbsArg.cxx:1817
RooAbsArg::treeNodeServerList
void treeNodeServerList(RooAbsCollection *list, const RooAbsArg *arg=0, Bool_t doBranch=kTRUE, Bool_t doLeaf=kTRUE, Bool_t valueOnly=kFALSE, Bool_t recurseNonDerived=kFALSE) const
Fill supplied list with nodes of the arg tree, following all server links, starting with ourself as t...
Definition: RooAbsArg.cxx:536
x
Double_t x[n]
Definition: legend1.C:17
RooAbsReal::_forceNumInt
Bool_t _forceNumInt
Definition: RooAbsReal.h:480
RooAbsPdf::getGenerator
virtual Int_t getGenerator(const RooArgSet &directVars, RooArgSet &generateVars, Bool_t staticInitOK=kTRUE) const
Load generatedVars with the subset of directVars that we can generate events for, and return a code t...
Definition: RooAbsPdf.cxx:2474
TString::Length
Ssiz_t Length() const
Definition: TString.h:410
coutI
#define coutI(a)
Definition: RooMsgService.h:30
indent
static void indent(ostringstream &buf, int indent_level)
Definition: TClingCallFunc.cxx:87
RooAbsReal
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:61
RooArgList::at
RooAbsArg * at(Int_t idx) const
Return object at given index, or nullptr if index is out of range.
Definition: RooArgList.h:72
RooLinkedList::fwdIterator
RooFIter fwdIterator() const
Create a one-time-use forward iterator for this list.
Definition: RooLinkedList.cxx:750
RooAbsPdf::initGenerator
virtual void initGenerator(Int_t code)
Interface for one-time initialization to setup the generator for the specified code.
Definition: RooAbsPdf.cxx:2484
RooProdPdf::groupProductTerms
void groupProductTerms(std::list< std::vector< RooArgSet * >> &groupedTerms, RooArgSet &outerIntDeps, const RooLinkedList &terms, const RooLinkedList &norms, const RooLinkedList &imps, const RooLinkedList &ints, const RooLinkedList &cross) const
Group product into terms that can be calculated independently.
Definition: RooProdPdf.cxx:1459
RooCacheManager::selectFromSet2
RooArgSet selectFromSet2(RooArgSet const &argSet, int index) const
Create RooArgSet contatining the objects that are both in the cached set 2.
Definition: RooCacheManager.h:333
RooAbsArg::setExpensiveObjectCache
virtual void setExpensiveObjectCache(RooExpensiveObjectCache &cache)
Definition: RooAbsArg.h:503
RooAbsCacheElement
RooAbsCacheElement is the abstract base class for objects to be stored in RooAbsCache cache manager o...
Definition: RooAbsCacheElement.h:26
RooProdPdf::_defNormSet
RooArgSet _defNormSet
Definition: RooProdPdf.h:182
TString
Basic string class.
Definition: TString.h:136
RooAbsCollection::GetName
const char * GetName() const
Returns name of object.
Definition: RooAbsCollection.h:286
RooAbsArg::expensiveObjectCache
RooExpensiveObjectCache & expensiveObjectCache() const
Definition: RooAbsArg.cxx:2277
RooGenProdProj
RooGenProdProj is an auxiliary class for RooProdPdf that calculates a general normalised projection o...
Definition: RooGenProdProj.h:26
RooCustomizer::build
RooAbsArg * build(const char *masterCatState, Bool_t verbose=kFALSE)
Build a clone of the prototype executing all registered 'replace' rules and 'split' rules for the mas...
Definition: RooCustomizer.cxx:397
TObject::InheritsFrom
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:445
v
@ v
Definition: rootcling_impl.cxx:3664
RooProdPdf::redirectServersHook
Bool_t redirectServersHook(const RooAbsCollection &, Bool_t, Bool_t, Bool_t)
Implement support for node removal.
Definition: RooProdPdf.cxx:2318
RooAbsPdf::isDirectGenSafe
virtual Bool_t isDirectGenSafe(const RooAbsArg &arg) const
Check if given observable can be safely generated using the pdfs internal generator mechanism (if tha...
Definition: RooAbsPdf.cxx:2509
RooProdPdf::initializeFromCmdArgList
void initializeFromCmdArgList(const RooArgSet &fullPdfSet, const RooLinkedList &l)
Initialize RooProdPdf configuration from given list of RooCmdArg configuration arguments and set of '...
Definition: RooProdPdf.cxx:374
RooListProxy::remove
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE) override
Reimplementation of standard RooArgList::remove()
Definition: RooListProxy.cxx:143
RooProdPdf::evaluateSpan
virtual RooSpan< double > evaluateSpan(RooBatchCompute::RunContext &evalData, const RooArgSet *normSet) const
Evaluate product of PDFs using input data in evalData.
Definition: RooProdPdf.cxx:509
bool
RooAbsRealLValue::getMin
virtual Double_t getMin(const char *name=0) const
Get miniminum of currently defined range.
Definition: RooAbsRealLValue.h:86
RooCustomizer::replaceArg
void replaceArg(const RooAbsArg &orig, const RooAbsArg &subst)
Replace any occurence of arg 'orig' with arg 'subst'.
Definition: RooCustomizer.cxx:341
RooAbsCollection::iterator
RooLinkedListIter iterator(Bool_t dir=kIterForward) const
TIterator-style iteration over contained elements.
Definition: RooAbsCollection.h:186
RooSpan::size
constexpr std::span< T >::index_type size() const noexcept
Definition: RooSpan.h:121
RooRangeBoolean.h
RooRealIntegral
RooRealIntegral performs hybrid numerical/analytical integrals of RooAbsReal objects.
Definition: RooRealIntegral.h:34
RooProdPdf::extendMode
virtual ExtendMode extendMode() const
If this product contains exactly one extendable p.d.f return the extension abilities of that p....
Definition: RooProdPdf.cxx:1852
RooGenProdProj.h
RooAbsCollection::setName
void setName(const char *name)
Definition: RooAbsCollection.h:282
RooFit::LinkStateMgmt
@ LinkStateMgmt
Definition: RooGlobalFunc.h:60
RooArgSet::snapshot
RooArgSet * snapshot(bool deepCopy=true) const
Use RooAbsCollection::snapshot(), but return as RooArgSet.
Definition: RooArgSet.h:115
RooProdPdf::CacheElem::_numList
RooArgList _numList
Definition: RooProdPdf.h:143
RooLinkedList::At
TObject * At(int index) const
Return object stored in sequential position given by index.
Definition: RooLinkedList.cxx:487
RooProdPdf::CacheElem
Definition: RooProdPdf.h:137
RooAbsPdf::ExtendMode
ExtendMode
Definition: RooAbsPdf.h:230
RooProdPdf::getCacheElem
CacheElem * getCacheElem(RooArgSet const *nset) const
Definition: RooProdPdf.cxx:448
RooTrace.h
group
Definition: TWinNTSystem.h:59
RooProdPdf::specializeRatio
RooAbsReal * specializeRatio(RooFormulaVar &input, const char *targetRangeName) const
Definition: RooProdPdf.cxx:1413
RooAICRegistry::store
Int_t store(const std::vector< Int_t > &codeList, RooArgSet *set1=0, RooArgSet *set2=0, RooArgSet *set3=0, RooArgSet *set4=0)
Store given arrays of integer codes, and up to four RooArgSets in the registry (each setX pointer may...
Definition: RooAICRegistry.cxx:100
RooAbsArg::setStringAttribute
void setStringAttribute(const Text_t *key, const Text_t *value)
Associate string 'value' to this object under key 'key'.
Definition: RooAbsArg.cxx:323
RooProduct
A RooProduct represents the product of a given set of RooAbsReal objects.
Definition: RooProduct.h:29
RooProdPdf::CacheElem::_isRearranged
Bool_t _isRearranged
Definition: RooProdPdf.h:147
RooProdPdf::CacheElem::_partList
RooArgList _partList
Definition: RooProdPdf.h:142
RooHelpers::getColonSeparatedNameString
std::string getColonSeparatedNameString(RooArgSet const &argSet)
Create a string with all sorted names of RooArgSet elements separated by colons.
Definition: RooHelpers.cxx:271
RooAbsCollection::contains
Bool_t contains(const RooAbsArg &var) const
Check if collection contains an argument with the same name as var.
Definition: RooAbsCollection.h:143
RooProdPdf::RooProdGenContext
friend class RooProdGenContext
Definition: RooProdPdf.h:164
RooSpan::end
constexpr std::span< T >::iterator end() const
Definition: RooSpan.h:102
RooCacheManager::getObjByIndex
T * getObjByIndex(Int_t index) const
Retrieve payload object by slot index.
Definition: RooCacheManager.h:310
RooLinkedListIter
A wrapper around TIterator derivatives.
Definition: RooLinkedListIter.h:204
RooFormulaVar
A RooFormulaVar is a generic implementation of a real-valued object, which takes a RooArgList of serv...
Definition: RooFormulaVar.h:30
RooAbsPdf::CanNotBeExtended
@ CanNotBeExtended
Definition: RooAbsPdf.h:230
RooCustomizer.h
RooProduct::components
RooArgList components()
Definition: RooProduct.h:47
RooProdPdf.h
TRACE_DESTROY
#define TRACE_DESTROY
Definition: RooTrace.h:24
RooProdPdf::CacheElem::_ownedList
RooArgList _ownedList
Definition: RooProdPdf.h:145
RooProdPdf::makeRGPPName
std::string makeRGPPName(const char *pfx, const RooArgSet &term, const RooArgSet &iset, const RooArgSet &nset, const char *isetRangeName) const
Make an appropriate automatic name for a RooGenProdProj object in getPartIntList()
Definition: RooProdPdf.cxx:1736
RooFormulaVar::getParameter
RooAbsArg * getParameter(const char *name) const
Return pointer to parameter with given name.
Definition: RooFormulaVar.h:42
RooProdPdf::_pdfNSetList
RooLinkedList _pdfNSetList
Definition: RooProdPdf.h:173
RooSpan::empty
constexpr bool empty() const noexcept
Definition: RooSpan.h:125
RooAbsReal::integralNameSuffix
TString integralNameSuffix(const RooArgSet &iset, const RooArgSet *nset=0, const char *rangeName=0, Bool_t omitEmpty=kFALSE) const
Construct string with unique suffix name to give to integral object that encodes integrated observabl...
Definition: RooAbsReal.cxx:792
RooProdPdf::initGenerator
virtual void initGenerator(Int_t code)
Forward one-time initialization call to component generation initialization methods.
Definition: RooProdPdf.cxx:1933
TNamed
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
RooProdPdf::_extendedIndex
Int_t _extendedIndex
Definition: RooProdPdf.h:174
RooFit::Generation
@ Generation
Definition: RooGlobalFunc.h:60
RooFIter
A one-time forward iterator working on RooLinkedList or RooAbsCollection.
Definition: RooLinkedListIter.h:40
RooProdPdf::setCacheAndTrackHints
virtual void setCacheAndTrackHints(RooArgSet &)
Label OK'ed components of a RooProdPdf with cache-and-track, and label all RooProdPdf descendants wit...
Definition: RooProdPdf.cxx:2252
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:101
RooAbsCollection::get
Storage_t const & get() const
Const access to the underlying stl container.
Definition: RooAbsCollection.h:76
RooProdPdf::findPdfNSet
RooArgSet * findPdfNSet(RooAbsPdf &pdf) const
Look up user specified normalization set for given input PDF component.
Definition: RooProdPdf.cxx:2044
RooLinkedList
RooLinkedList is an collection class for internal use, storing a collection of RooAbsArg pointers in ...
Definition: RooLinkedList.h:37
RooProdPdf::CacheElem::_rearrangedDen
std::unique_ptr< RooAbsReal > _rearrangedDen
Definition: RooProdPdf.h:149
RooAbsPdf::normRange
const char * normRange() const
Definition: RooAbsPdf.h:264
RooFIter::next
RooAbsArg * next()
Return next element or nullptr if at end.
Definition: RooLinkedListIter.h:49
RooRealProxy.h
RooAddition
RooAddition calculates the sum of a set of RooAbsReal terms, or when constructed with two sets,...
Definition: RooAddition.h:27
RooAbsCollection
RooAbsCollection is an abstract container object that can hold multiple RooAbsArg objects.
Definition: RooAbsCollection.h:33
RooAbsArg::RooArgSet
friend class RooArgSet
Definition: RooAbsArg.h:599
RooAbsGenContext
RooAbsGenContext is the abstract base class for generator contexts of RooAbsPdf objects.
Definition: RooAbsGenContext.h:26
TRACE_CREATE
#define TRACE_CREATE
Definition: RooTrace.h:23
RooAbsCollection::add
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
Definition: RooAbsCollection.cxx:455
RooAddition.h
RooAbsArg::getStringAttribute
const Text_t * getStringAttribute(const Text_t *key) const
Get string attribute mapped under key 'key'.
Definition: RooAbsArg.cxx:336
RooAbsCollection::size
Storage_t::size_type size() const
Definition: RooAbsCollection.h:214
RooProdPdf::genContext
virtual RooAbsGenContext * genContext(const RooArgSet &vars, const RooDataSet *prototype=0, const RooArgSet *auxProto=0, Bool_t verbose=kFALSE) const
Return generator context optimized for generating events from product p.d.f.s.
Definition: RooProdPdf.cxx:1878
RooLinkedList::Add
virtual void Add(TObject *arg)
Definition: RooLinkedList.h:64
RooAbsCollection::printValue
virtual void printValue(std::ostream &os) const
Print value of collection, i.e.
Definition: RooAbsCollection.cxx:1093
RooLinkedList::iterator
RooLinkedListIter iterator(Bool_t forward=kTRUE) const
Create an iterator for this list.
Definition: RooLinkedList.cxx:741
RooProdPdf::clone
virtual TObject * clone(const char *newname) const
Definition: RooProdPdf.h:60
RooSpan::begin
constexpr std::span< T >::iterator begin() const
Definition: RooSpan.h:98
RooRealVar.h
RooProdPdf::getAnalyticalIntegralWN
Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &numVars, const RooArgSet *normSet, const char *rangeName=0) const
Determine which part (if any) of given integral can be performed analytically.
Definition: RooProdPdf.cxx:1785
RooAbsCollection::overlaps
bool overlaps(Iterator_t otherCollBegin, Iterator_t otherCollEnd) const
Definition: RooAbsCollection.h:160
RooProduct.h
RooAbsCollection::addOwned
virtual Bool_t addOwned(RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
Definition: RooAbsCollection.cxx:403
RooRealIntegral::intRange
const char * intRange()
Definition: RooRealIntegral.h:57
RooProdPdf::rearrangeProduct
void rearrangeProduct(CacheElem &) const
Definition: RooProdPdf.cxx:1143
RooProdPdf::analyticalIntegralWN
Double_t analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=0) const
Return analytical integral defined by given scenario code.
Definition: RooProdPdf.cxx:1805
RooConstVar.h
RooAbsCollection::selectCommon
bool selectCommon(const RooAbsCollection &refColl, RooAbsCollection &outColl) const
Create a subset of the current collection, consisting only of those elements that are contained as we...
Definition: RooAbsCollection.cxx:702
RooProdPdf::processProductTerm
std::vector< RooAbsReal * > processProductTerm(const RooArgSet *nset, const RooArgSet *iset, const char *isetRangeName, const RooArgSet *term, const RooArgSet &termNSet, const RooArgSet &termISet, Bool_t &isOwned, Bool_t forceWrap=kFALSE) const
Calculate integrals of factorized product terms over observables iset while normalized to observables...
Definition: RooProdPdf.cxx:1550
RooGlobalFunc.h
RooProdPdf::expectedEvents
virtual Double_t expectedEvents(const RooArgSet *nset) const
Return the expected number of events associated with the extendable input PDF in the product.
Definition: RooProdPdf.cxx:1863
Double_t
double Double_t
Definition: RtypesCore.h:59
RooAbsArg::getObservables
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Given a set of possible observables, return the observables that this PDF depends on.
Definition: RooAbsArg.h:295
RooAbsArg::addOwnedComponents
Bool_t addOwnedComponents(const RooArgSet &comps)
Take ownership of the contents of 'comps'.
Definition: RooAbsArg.cxx:2216
RooProdPdf::CacheElem::_denList
RooArgList _denList
Definition: RooProdPdf.h:144
RooProdPdf::binBoundaries
virtual std::list< Double_t > * binBoundaries(RooAbsRealLValue &, Double_t, Double_t) const
Forward the plot sampling hint from the p.d.f. that defines the observable obs.
Definition: RooProdPdf.cxx:2232
RooProdPdf::RooProdPdf
RooProdPdf()
Default constructor.
Definition: RooProdPdf.cxx:81
RooAbsCacheElement::Action
Action
Definition: RooAbsCacheElement.h:39
RooProdPdf::~RooProdPdf
virtual ~RooProdPdf()
Destructor.
Definition: RooProdPdf.cxx:441
TObject
Mother of all ROOT objects.
Definition: TObject.h:37
coutF
#define coutF(a)
Definition: RooMsgService.h:34
RooAbsArg::setAttribute
void setAttribute(const Text_t *name, Bool_t value=kTRUE)
Set (default) or clear a named boolean attribute of this object.
Definition: RooAbsArg.cxx:291
RooProdPdf::CacheElem::containedArgs
RooArgList containedArgs(Action) override
Return RooAbsArg components contained in the cache.
Definition: RooProdPdf.cxx:1977
RooNameReg::str
static const char * str(const TNamed *ptr)
Return C++ string corresponding to given TNamed pointer.
Definition: RooNameReg.cxx:103
RooProdPdf::checkObservables
virtual Bool_t checkObservables(const RooArgSet *nset) const
Obsolete.
Definition: RooProdPdf.cxx:1842
name
char name[80]
Definition: TGX11.cxx:110
RooAbsArg::dependsOn
Bool_t dependsOn(const RooAbsCollection &serverList, const RooAbsArg *ignoreArg=0, Bool_t valueOnly=kFALSE) const
Test whether we depend on (ie, are served by) any object in the specified collection.
Definition: RooAbsArg.cxx:799
RooProdPdf::calculate
Double_t calculate(const RooProdPdf::CacheElem &cache, Bool_t verbose=kFALSE) const
Calculate running product of pdfs terms, using the supplied normalization set in 'normSetList' for ea...
Definition: RooProdPdf.cxx:475
ROOT::Math::Chebyshev::T
double T(double x)
Definition: ChebyshevPol.h:34
RooProdPdf::_useDefaultGen
Bool_t _useDefaultGen
Definition: RooProdPdf.h:177
RooProdPdf::getGenerator
virtual Int_t getGenerator(const RooArgSet &directVars, RooArgSet &generateVars, Bool_t staticInitOK=kTRUE) const
Query internal generation capabilities of component p.d.f.s and aggregate capabilities into master co...
Definition: RooProdPdf.cxx:1891
genreflex::verbose
bool verbose
Definition: rootcling_impl.cxx:133
RooAddition::list1
const RooArgList & list1() const
Definition: RooAddition.h:42
RooDataSet
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:33
RooCacheManager::lastIndex
Int_t lastIndex() const
Definition: RooCacheManager.h:66
dologD
#define dologD(a)
Definition: RooMsgService.h:65
RooAbsArg
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:72
RooAbsReal::getValues
virtual RooSpan< const double > getValues(RooBatchCompute::RunContext &evalData, const RooArgSet *normSet=nullptr) const
by this change, please consult the release notes for ROOT 6.24 for guidance on how to make this trans...
Definition: RooAbsReal.cxx:312
RooListProxy::add
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE) override
Reimplementation of standard RooArgList::add()
Definition: RooListProxy.cxx:104
RooRealIntegral.h
RooAbsPdf
Definition: RooAbsPdf.h:41
TNamed::GetName
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
RooProdPdf::generateEvent
virtual void generateEvent(Int_t code)
Generate a single event with configuration specified by 'code' Defer internal generation to component...
Definition: RooProdPdf.cxx:1956
RooProdPdf::fixRefRange
void fixRefRange(const char *rangeName)
Definition: RooProdPdf.cxx:2182
RooHelpers
Definition: RooHelpers.h:33
RooProdPdf::printMetaArgs
void printMetaArgs(std::ostream &os) const
Customized printing of arguments of a RooProdPdf to more intuitively reflect the contents of the prod...
Definition: RooProdPdf.cxx:2284
RooAICRegistry::retrieve
const std::vector< Int_t > & retrieve(Int_t masterCode) const
Retrieve the array of integer codes associated with the given master code.
Definition: RooAICRegistry.cxx:153
RooProdPdf::makeCondPdfRatioCorr
RooAbsReal * makeCondPdfRatioCorr(RooAbsReal &term, const RooArgSet &termNset, const RooArgSet &termImpSet, const char *normRange, const char *refRange) const
For single normalization ranges.
Definition: RooProdPdf.cxx:1126
Class
void Class()
Definition: Class.C:29
RooRealVar
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:39
RooAbsCollection::removeAll
virtual void removeAll()
Remove all arguments from our set, deleting them if we own them.
Definition: RooAbsCollection.cxx:644
RooFit::Eval
@ Eval
Definition: RooGlobalFunc.h:61
RooRealIntegral::integrand
const RooAbsReal & integrand() const
Definition: RooRealIntegral.h:58
RooAbsRealLValue
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
Definition: RooAbsRealLValue.h:31
RooProdPdf
RooProdPdf is an efficient implementation of a product of PDFs of the form.
Definition: RooProdPdf.h:37
RooRealIntegral::intVars
RooArgSet intVars() const
Definition: RooRealIntegral.h:56
RooProdPdf::CacheElem::writeToStream
void writeToStream(std::ostream &os) const
Definition: RooProdPdf.cxx:2337
RooLinkedList::IndexOf
Int_t IndexOf(const char *name) const
Return position of given object in list.
Definition: RooLinkedList.cxx:700
RooNameReg.h
RooAbsPdf::RooRealIntegral
friend class RooRealIntegral
Definition: RooAbsPdf.h:317
RooBatchCompute::RunContext
This struct enables passing computation data around between elements of a computation graph.
Definition: RunContext.h:31
RooCacheManager::selectFromSet1
RooArgSet selectFromSet1(RooArgSet const &argSet, int index) const
Create RooArgSet contatining the objects that are both in the cached set 1.
Definition: RooCacheManager.h:324
RooAbsPdf::canBeExtended
Bool_t canBeExtended() const
If true, PDF can provide extended likelihood term.
Definition: RooAbsPdf.h:236
RooLinkedList::Delete
void Delete(Option_t *o=0)
Remove all elements in collection and delete all elements NB: Collection does not own elements,...
Definition: RooLinkedList.cxx:575
RooAbsCollection::getSize
Int_t getSize() const
Definition: RooAbsCollection.h:231
RooSpan
A simple container to hold a batch of data values.
Definition: RooSpan.h:34
RooArgSet
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:33
RooAbsCollection::index
Int_t index(const RooAbsArg *arg) const
Returns index of given arg, or -1 if arg is not in the collection.
Definition: RooAbsCollection.h:247
RooLinkedList::Remove
virtual Bool_t Remove(TObject *arg)
Remove object from collection.
Definition: RooLinkedList.cxx:446
cxcoutD
#define cxcoutD(a)
Definition: RooMsgService.h:81
RooProdPdf::getParametersHook
virtual void getParametersHook(const RooArgSet *, RooArgSet *, Bool_t stripDisconnected) const
Definition: RooProdPdf.cxx:2133
int
RooProdPdf::factorizeProduct
void factorizeProduct(const RooArgSet &normSet, const RooArgSet &intSet, RooLinkedList &termList, RooLinkedList &normList, RooLinkedList &impDepList, RooLinkedList &crossDepList, RooLinkedList &intList) const
Factorize product in irreducible terms for given choice of integration/normalization.
Definition: RooProdPdf.cxx:598
RooAbsArg::getAttribute
Bool_t getAttribute(const Text_t *name) const
Check if a named attribute is set. By default, all attributes are unset.
Definition: RooAbsArg.cxx:314
RooFit::RooConst
RooConstVar & RooConst(Double_t val)
Definition: RooGlobalFunc.cxx:383