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