Logo ROOT   6.08/07
Reference Guide
RooRealIntegral.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 RooRealIntegral.cxx
19 \class RooRealIntegral
20 \ingroup Roofitcore
21 
22 RooRealIntegral performs hybrid numerical/analytical integrals of RooAbsReal objects
23 The class performs none of the actual integration, but only manages the logic
24 of what variables can be integrated analytically, accounts for eventual jacobian
25 terms and defines what numerical integrations needs to be done to complement the
26 analytical integral.
27 The actual analytical integrations (if any) are done in the PDF themselves, the numerical
28 integration is performed in the various implemenations of the RooAbsIntegrator base class.
29 **/
30 
31 #include "RooFit.h"
32 
33 #include "TClass.h"
34 #include "RooMsgService.h"
35 #include "Riostream.h"
36 #include "TObjString.h"
37 #include "TH1.h"
38 #include "RooRealIntegral.h"
39 #include "RooArgSet.h"
40 #include "RooAbsRealLValue.h"
41 #include "RooAbsCategoryLValue.h"
42 #include "RooRealBinding.h"
43 #include "RooRealAnalytic.h"
44 #include "RooInvTransform.h"
45 #include "RooSuperCategory.h"
46 #include "RooNumIntFactory.h"
47 #include "RooNumIntConfig.h"
48 #include "RooNameReg.h"
50 #include "RooConstVar.h"
51 #include "RooDouble.h"
52 #include "RooTrace.h"
53 
54 using namespace std;
55 
57 ;
58 
59 
61 
62 
63 ////////////////////////////////////////////////////////////////////////////////
64 
66  _valid(kFALSE),
67  _funcNormSet(0),
68  _iconfig(0),
69  _sumCatIter(0),
70  _mode(0),
71  _intOperMode(Hybrid),
72  _restartNumIntEngine(kFALSE),
73  _numIntEngine(0),
74  _numIntegrand(0),
75  _rangeName(0),
76  _params(0),
77  _cacheNum(kFALSE)
78 {
82 }
83 
84 
85 
86 ////////////////////////////////////////////////////////////////////////////////
87 /// Construct integral of 'function' over observables in 'depList'
88 /// in range 'rangeName' with normalization observables 'funcNormSet'
89 /// (for p.d.f.s). In the integral is performed to the maximum extent
90 /// possible the internal (analytical) integrals advertised by function.
91 /// The other integrations are performed numerically. The optional
92 /// config object prescribes how these numeric integrations are configured.
93 ///
94 
95 RooRealIntegral::RooRealIntegral(const char *name, const char *title,
96  const RooAbsReal& function, const RooArgSet& depList,
97  const RooArgSet* funcNormSet, const RooNumIntConfig* config,
98  const char* rangeName) :
99  RooAbsReal(name,title),
100  _valid(kTRUE),
101  _sumList("!sumList","Categories to be summed numerically",this,kFALSE,kFALSE),
102  _intList("!intList","Variables to be integrated numerically",this,kFALSE,kFALSE),
103  _anaList("!anaList","Variables to be integrated analytically",this,kFALSE,kFALSE),
104  _jacList("!jacList","Jacobian product term",this,kFALSE,kFALSE),
105  _facList("!facList","Variables independent of function",this,kFALSE,kTRUE),
106  _facListIter(_facList.createIterator()),
107  _jacListIter(_jacList.createIterator()),
108  _function("!func","Function to be integrated",this,
109  const_cast<RooAbsReal&>(function),kFALSE,kFALSE),
110  _iconfig((RooNumIntConfig*)config),
111  _sumCat("!sumCat","SuperCategory for summation",this,kFALSE,kFALSE),
112  _sumCatIter(0),
113  _mode(0),
116  _numIntEngine(0),
117  _numIntegrand(0),
118  _rangeName((TNamed*)RooNameReg::ptr(rangeName)),
119  _params(0),
121 {
122  // A) Check that all dependents are lvalues
123  //
124  // B) Check if list of dependents can be re-expressed in
125  // lvalues that are higher in the expression tree
126  //
127  // C) Check for dependents that the PDF insists on integrating
128  // analytically iself
129  //
130  // D) Make list of servers that can be integrated analytically
131  // Add all parameters/dependents as value/shape servers
132  //
133  // E) Interact with function to make list of objects actually integrated analytically
134  //
135  // F) Make list of numerical integration variables consisting of:
136  // - Category dependents of RealLValues in analytical integration
137  // - Leaf nodes server lists of function server that are not analytically integrated
138  // - Make Jacobian list for analytically integrated RealLValues
139  //
140  // G) Split numeric list in integration list and summation list
141  //
142 
143  oocxcoutI(&function,Integration) << "RooRealIntegral::ctor(" << GetName() << ") Constructing integral of function "
144  << function.GetName() << " over observables" << depList << " with normalization "
145  << (funcNormSet?*funcNormSet:RooArgSet()) << " with range identifier "
146  << (rangeName?rangeName:"<none>") << endl ;
147 
148 
149  // Choose same expensive object cache as integrand
151 // cout << "RRI::ctor(" << GetName() << ") setting expensive object cache to " << &expensiveObjectCache() << " as taken from " << function.GetName() << endl ;
152 
153  // Use objects integrator configuration if none is specified
154  if (!_iconfig) _iconfig = (RooNumIntConfig*) function.getIntegratorConfig() ;
155 
156  // Save private copy of funcNormSet, if supplied, excluding factorizing terms
157  if (funcNormSet) {
158  _funcNormSet = new RooArgSet ;
159  TIterator* iter = funcNormSet->createIterator() ;
160  RooAbsArg* nArg ;
161  while ((nArg=(RooAbsArg*)iter->Next())) {
162  if (function.dependsOn(*nArg)) {
163  _funcNormSet->addClone(*nArg) ;
164  }
165  }
166  delete iter ;
167  } else {
168  _funcNormSet = 0 ;
169  }
170 
171  //_funcNormSet = funcNormSet ? (RooArgSet*)funcNormSet->snapshot(kFALSE) : 0 ;
172 
173  // Make internal copy of dependent list
174  RooArgSet intDepList(depList) ;
175 
176  RooAbsArg *arg ;
177 
178  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
179  // * A) Check that all dependents are lvalues and filter out any
180  // dependents that the PDF doesn't explicitly depend on
181  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
182 
183  TIterator* depIter = intDepList.createIterator() ;
184  while((arg=(RooAbsArg*)depIter->Next())) {
185  if(!arg->isLValue()) {
186  coutE(InputArguments) << ClassName() << "::" << GetName() << ": cannot integrate non-lvalue ";
187  arg->Print("1");
188  _valid= kFALSE;
189  }
190  if (!function.dependsOn(*arg)) {
191  RooAbsArg* argClone = (RooAbsArg*) arg->Clone() ;
192  _facListOwned.addOwned(*argClone) ;
193  _facList.add(*argClone) ;
194  addServer(*argClone,kFALSE,kTRUE) ;
195  }
196  }
197 
198  if (_facList.getSize()>0) {
199  oocxcoutI(&function,Integration) << function.GetName() << ": Factorizing obserables are " << _facList << endl ;
200  }
201 
202 
203  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
204  // * B) Check if list of dependents can be re-expressed in *
205  // * lvalues that are higher in the expression tree *
206  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
207 
208 
209  // Initial fill of list of LValue branches
210  RooArgSet exclLVBranches("exclLVBranches") ;
211  RooArgSet branchList,branchListVD ;
212  function.branchNodeServerList(&branchList) ;
213 
214  TIterator* bIter = branchList.createIterator() ;
215  RooAbsArg* branch ;
216  while((branch=(RooAbsArg*)bIter->Next())) {
217  RooAbsRealLValue *realArgLV = dynamic_cast<RooAbsRealLValue*>(branch) ;
218  RooAbsCategoryLValue *catArgLV = dynamic_cast<RooAbsCategoryLValue*>(branch) ;
219  if ((realArgLV && (realArgLV->isJacobianOK(intDepList)!=0)) || catArgLV) {
220  exclLVBranches.add(*branch) ;
221 // cout << "exclv branch = " << endl ;
222 // branch->printCompactTree() ;
223  }
224  if (dependsOnValue(*branch)) {
225  branchListVD.add(*branch) ;
226  } else {
227 // cout << "value of self does not depend on branch " << branch->GetName() << endl ;
228  }
229  }
230  delete bIter ;
231  exclLVBranches.remove(depList,kTRUE,kTRUE) ;
232 // cout << "exclLVBranches = " << exclLVBranches << endl ;
233 
234  // Initial fill of list of LValue leaf servers (put in intDepList)
235  RooArgSet exclLVServers("exclLVServers") ;
236  exclLVServers.add(intDepList) ;
237 
238 // cout << "begin exclLVServers = " << exclLVServers << endl ;
239 
240  // Obtain mutual exclusive dependence by iterative reduction
241  TIterator *sIter = exclLVServers.createIterator() ;
242  bIter = exclLVBranches.createIterator() ;
243  RooAbsArg *server ;
244  Bool_t converged(kFALSE) ;
245  while(!converged) {
246  converged=kTRUE ;
247 
248  // Reduce exclLVServers to only those serving exclusively exclLVBranches
249  sIter->Reset() ;
250  while ((server=(RooAbsArg*)sIter->Next())) {
251  if (!servesExclusively(server,exclLVBranches,branchListVD)) {
252  exclLVServers.remove(*server) ;
253 // cout << "removing " << server->GetName() << " from exclLVServers because servesExclusively(" << server->GetName() << "," << exclLVBranches << ") faile" << endl ;
254  converged=kFALSE ;
255  }
256  }
257 
258  // Reduce exclLVBranches to only those depending exclusisvely on exclLVservers
259  bIter->Reset() ;
260  while((branch=(RooAbsArg*)bIter->Next())) {
261  RooArgSet* brDepList = branch->getObservables(&intDepList) ;
262  RooArgSet bsList(*brDepList,"bsList") ;
263  delete brDepList ;
264  bsList.remove(exclLVServers,kTRUE,kTRUE) ;
265  if (bsList.getSize()>0) {
266  exclLVBranches.remove(*branch,kTRUE,kTRUE) ;
267 // cout << "removing " << branch->GetName() << " from exclLVBranches" << endl ;
268  converged=kFALSE ;
269  }
270  }
271  }
272 
273  // Eliminate exclLVBranches that do not depend on any LVServer
274  bIter->Reset() ;
275  while((branch=(RooAbsArg*)bIter->Next())) {
276  if (!branch->dependsOnValue(exclLVServers)) {
277  //cout << "LV branch " << branch->GetName() << " does not depend on any LVServer (" << exclLVServers << ") and will be removed" << endl ;
278  exclLVBranches.remove(*branch,kTRUE,kTRUE) ;
279  }
280  }
281 
282  delete sIter ;
283  delete bIter ;
284 
285 // cout << "end exclLVServers = " << exclLVServers << endl ;
286 
287  // Replace exclusive lvalue branch servers with lvalue branches
288  // WVE Don't do this for binned distributions - deal with this using numeric integration with transformed bin boundaroes
289  if (exclLVServers.getSize()>0 && !function.isBinnedDistribution(exclLVBranches)) {
290 // cout << "activating LVservers " << exclLVServers << " for use in integration " << endl ;
291  intDepList.remove(exclLVServers) ;
292  intDepList.add(exclLVBranches) ;
293 
294  //cout << "intDepList removing exclLVServers " << exclLVServers << endl ;
295  //cout << "intDepList adding exclLVBranches " << exclLVBranches << endl ;
296 
297  }
298 
299 
300  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
301  // * C) Check for dependents that the PDF insists on integrating *
302  // analytically iself *
303  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
304 
305  RooArgSet anIntOKDepList ;
306  depIter->Reset() ;
307  while((arg=(RooAbsArg*)depIter->Next())) {
308  if (function.forceAnalyticalInt(*arg)) {
309  anIntOKDepList.add(*arg) ;
310  }
311  }
312  delete depIter ;
313 
314  if (anIntOKDepList.getSize()>0) {
315  oocxcoutI(&function,Integration) << function.GetName() << ": Observables that function forcibly requires to be integrated internally " << anIntOKDepList << endl ;
316  }
317 
318 
319  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
320  // * D) Make list of servers that can be integrated analytically *
321  // Add all parameters/dependents as value/shape servers *
322  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
323 
324  sIter = function.serverIterator() ;
325  while((arg=(RooAbsArg*)sIter->Next())) {
326 
327  //cout << "considering server" << arg->GetName() << endl ;
328 
329  // Dependent or parameter?
330  if (!arg->dependsOnValue(intDepList)) {
331 
332  //cout << " server does not depend on observables, adding server as value server to integral" << endl ;
333 
334  if (function.dependsOnValue(*arg)) {
335  addServer(*arg,kTRUE,kFALSE) ;
336  }
337 
338  continue ;
339 
340  } else {
341 
342  // Add final dependents of arg as shape servers
343  RooArgSet argLeafServers ;
344  arg->leafNodeServerList(&argLeafServers,0,kFALSE) ;
345 
346  //arg->printCompactTree() ;
347  //cout << "leaf nodes of server are " << argLeafServers << " depList = " << depList << endl ;
348 
349  // Skip arg if it is neither value or shape server
350  if (!arg->isValueServer(function) && !arg->isShapeServer(function)) {
351  //cout << " server is neither value not shape server of function, ignoring" << endl ;
352  continue ;
353  }
354 
355  TIterator* lIter = argLeafServers.createIterator() ;
356  RooAbsArg* leaf ;
357  while((leaf=(RooAbsArg*)lIter->Next())) {
358 
359  //cout << " considering leafnode " << leaf->GetName() << " of server " << arg->GetName() << endl ;
360 
361  if (depList.find(leaf->GetName()) && function.dependsOnValue(*leaf)) {
362 
363  RooAbsRealLValue* leaflv = dynamic_cast<RooAbsRealLValue*>(leaf) ;
364  if (leaflv && leaflv->getBinning(rangeName).isParameterized()) {
365  oocxcoutD(&function,Integration) << function.GetName() << " : Observable " << leaf->GetName() << " has parameterized binning, add value dependence of boundary objects rather than shape of leaf" << endl ;
366  if (leaflv->getBinning(rangeName).lowBoundFunc()) {
367  addServer(*leaflv->getBinning(rangeName).lowBoundFunc(),kTRUE,kFALSE) ;
368  }
369  if(leaflv->getBinning(rangeName).highBoundFunc()) {
370  addServer(*leaflv->getBinning(rangeName).highBoundFunc(),kTRUE,kFALSE) ;
371  }
372  } else {
373  oocxcoutD(&function,Integration) << function.GetName() << ": Adding observable " << leaf->GetName() << " of server "
374  << arg->GetName() << " as shape dependent" << endl ;
375  addServer(*leaf,kFALSE,kTRUE) ;
376  }
377  } else if (!depList.find(leaf->GetName())) {
378 
379  if (function.dependsOnValue(*leaf)) {
380  oocxcoutD(&function,Integration) << function.GetName() << ": Adding parameter " << leaf->GetName() << " of server " << arg->GetName() << " as value dependent" << endl ;
381  addServer(*leaf,kTRUE,kFALSE) ;
382  } else {
383  oocxcoutD(&function,Integration) << function.GetName() << ": Adding parameter " << leaf->GetName() << " of server " << arg->GetName() << " as shape dependent" << endl ;
384  addServer(*leaf,kFALSE,kTRUE) ;
385  }
386  }
387  }
388  delete lIter ;
389  }
390 
391  // If this dependent arg is self-normalized, stop here
392  //if (function.selfNormalized()) continue ;
393 
394  Bool_t depOK(kFALSE) ;
395  // Check for integratable AbsRealLValue
396 
397  //cout << "checking server " << arg->IsA()->GetName() << "::" << arg->GetName() << endl ;
398 
399  if (arg->isDerived()) {
400  RooAbsRealLValue *realArgLV = dynamic_cast<RooAbsRealLValue*>(arg) ;
401  RooAbsCategoryLValue *catArgLV = dynamic_cast<RooAbsCategoryLValue*>(arg) ;
402  //cout << "realArgLV = " << realArgLV << " intDepList = " << intDepList << endl ;
403  if ((realArgLV && intDepList.find(realArgLV->GetName()) && (realArgLV->isJacobianOK(intDepList)!=0)) || catArgLV) {
404 
405  //cout << " arg " << arg->GetName() << " is derived LValue with valid jacobian" << endl ;
406 
407  // Derived LValue with valid jacobian
408  depOK = kTRUE ;
409 
410  // Now, check for overlaps
411  Bool_t overlapOK = kTRUE ;
412  RooAbsArg *otherArg ;
413  TIterator* sIter2 = function.serverIterator() ;
414  while((otherArg=(RooAbsArg*)sIter2->Next())) {
415  // skip comparison with self
416  if (arg==otherArg) continue ;
417  if (otherArg->IsA()==RooConstVar::Class()) continue ;
418  if (arg->overlaps(*otherArg,kTRUE)) {
419  //cout << "arg " << arg->GetName() << " overlaps with " << otherArg->GetName() << endl ;
420  //overlapOK=kFALSE ;
421  }
422  }
423  // coverity[DEADCODE]
424  if (!overlapOK) depOK=kFALSE ;
425 
426  //cout << "overlap check returns OK=" << (depOK?"T":"F") << endl ;
427 
428  delete sIter2 ;
429  }
430  } else {
431  // Fundamental types are always OK
432  depOK = kTRUE ;
433  }
434 
435  // Add server to list of dependents that are OK for analytical integration
436  if (depOK) {
437  anIntOKDepList.add(*arg,kTRUE) ;
438  oocxcoutI(&function,Integration) << function.GetName() << ": Observable " << arg->GetName() << " is suitable for analytical integration (if supported by p.d.f)" << endl ;
439  }
440  }
441  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
442  // * E) interact with function to make list of objects actually integrated analytically *
443  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
444 
445  RooArgSet anIntDepList ;
446 
447  RooArgSet *anaSet = new RooArgSet( _anaList, Form("UniqueCloneOf_%s",_anaList.GetName()));
449  _anaList.removeAll() ;
450  _anaList.add(*anaSet);
451  delete anaSet;
452 
453  // Avoid confusion -- if mode is zero no analytical integral is defined regardless of contents of _anaListx
454  if (_mode==0) {
455  _anaList.removeAll() ;
456  }
457 
458  if (_mode!=0) {
459  oocxcoutI(&function,Integration) << function.GetName() << ": Function integrated observables " << _anaList << " internally with code " << _mode << endl ;
460  }
461 
462 
463  // WVE kludge: synchronize dset for use in analyticalIntegral
464  function.getVal(_funcNormSet) ;
465 
466  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
467  // * F) Make list of numerical integration variables consisting of: *
468  // * - Category dependents of RealLValues in analytical integration *
469  // * - Expanded server lists of server that are not analytically integrated *
470  // * Make Jacobian list with analytically integrated RealLValues *
471  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
472 
473  RooArgSet numIntDepList ;
474 
475 
476  // Loop over actually analytically integrated dependents
477  TIterator* aiIter = _anaList.createIterator() ;
478  while ((arg=(RooAbsArg*)aiIter->Next())) {
479 
480  // Process only derived RealLValues
481  if (arg->IsA()->InheritsFrom(RooAbsRealLValue::Class()) && arg->isDerived() && !arg->isFundamental()) {
482 
483  // Add to list of Jacobians to calculate
484  _jacList.add(*arg) ;
485 
486  // Add category dependent of LValueReal used in integration
487  RooAbsArg *argDep ;
488  RooArgSet *argDepList = arg->getObservables(&intDepList) ;
489  TIterator *adIter = argDepList->createIterator() ;
490  while ((argDep=(RooAbsArg*)adIter->Next())) {
491  if (argDep->IsA()->InheritsFrom(RooAbsCategoryLValue::Class()) && intDepList.contains(*argDep)) {
492  numIntDepList.add(*argDep,kTRUE) ;
493  }
494  }
495  delete adIter ;
496  delete argDepList ;
497  }
498  }
499  delete aiIter ;
500 
501 
502  // If nothing was integrated analytically, swap back LVbranches for LVservers for subsequent numeric integration
503  if (_anaList.getSize()==0) {
504  if (exclLVServers.getSize()>0) {
505  //cout << "NUMINT phase analList is empty. exclLVServers = " << exclLVServers << endl ;
506  intDepList.remove(exclLVBranches) ;
507  intDepList.add(exclLVServers) ;
508  }
509  }
510  //cout << "NUMINT intDepList = " << intDepList << endl ;
511 
512  // Loop again over function servers to add remaining numeric integrations
513  sIter->Reset() ;
514  while((arg=(RooAbsArg*)sIter->Next())) {
515 
516  //cout << "processing server for numeric integration " << arg->IsA()->GetName() << "::" << arg->GetName() << endl ;
517 
518  // Process only servers that are not treated analytically
519  if (!_anaList.find(arg->GetName()) && arg->dependsOn(intDepList)) {
520 
521  // Process only derived RealLValues
522  if (dynamic_cast<RooAbsLValue*>(arg) && arg->isDerived() && intDepList.contains(*arg)) {
523  numIntDepList.add(*arg,kTRUE) ;
524  } else {
525 
526  // WVE this will only get the observables, but not l-value transformations
527  // Expand server in final dependents
528  RooArgSet *argDeps = arg->getObservables(&intDepList) ;
529 
530  if (argDeps->getSize()>0) {
531 
532  // Add final dependents, that are not forcibly integrated analytically,
533  // to numerical integration list
534  TIterator* iter = argDeps->createIterator() ;
535  RooAbsArg* dep ;
536  while((dep=(RooAbsArg*)iter->Next())) {
537  if (!_anaList.find(dep->GetName())) {
538  numIntDepList.add(*dep,kTRUE) ;
539  }
540  }
541  delete iter ;
542 
543  }
544  delete argDeps ;
545  }
546 
547  }
548  }
549  delete sIter ;
550 
551  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
552  // * G) Split numeric list in integration list and summation list *
553  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
554 
555  // Split numeric integration list in summation and integration lists
556  TIterator* numIter=numIntDepList.createIterator() ;
557  while ((arg=(RooAbsArg*)numIter->Next())) {
558 
559  if (arg->IsA()->InheritsFrom(RooAbsRealLValue::Class())) {
560  _intList.add(*arg) ;
561  } else if (arg->IsA()->InheritsFrom(RooAbsCategoryLValue::Class())) {
562  _sumList.add(*arg) ;
563  }
564  }
565  delete numIter ;
566 
567  if (_anaList.getSize()>0) {
568  oocxcoutI(&function,Integration) << function.GetName() << ": Observables " << _anaList << " are analytically integrated with code " << _mode << endl ;
569  }
570  if (_intList.getSize()>0) {
571  oocxcoutI(&function,Integration) << function.GetName() << ": Observables " << _intList << " are numerically integrated" << endl ;
572  }
573  if (_sumList.getSize()>0) {
574  oocxcoutI(&function,Integration) << function.GetName() << ": Observables " << _sumList << " are numerically summed" << endl ;
575  }
576 
577 
578  // Determine operating mode
579  if (numIntDepList.getSize()>0) {
580  // Numerical and optional Analytical integration
581  _intOperMode = Hybrid ;
582  } else if (_anaList.getSize()>0) {
583  // Purely analytical integration
585  } else {
586  // No integration performed
588  }
589 
590  // Determine auto-dirty status
592 
593  // Create value caches for _intList and _sumList
596 
597 
598  if (_sumList.getSize()>0) {
599  RooSuperCategory *sumCat = new RooSuperCategory(Form("%s_sumCat",GetName()),"sumCat",_sumList) ;
600  _sumCatIter = sumCat->typeIterator() ;
601  _sumCat.addOwned(*sumCat) ;
602  }
603 
605 }
606 
607 
608 
609 ////////////////////////////////////////////////////////////////////////////////
610 /// Set appropriate cache operation mode for integral depending on cache operation
611 /// mode of server objects
612 
614 {
615  // If any of our servers are is forcedDirty or a projectedDependent, then we need to be ADirty
616  TIterator* siter = serverIterator() ;
617  RooAbsArg* server ;
618  while((server=(RooAbsArg*)siter->Next())){
619  if (server->isValueServer(*this)) {
620  RooArgSet leafSet ;
621  server->leafNodeServerList(&leafSet) ;
622  TIterator* liter = leafSet.createIterator() ;
623  RooAbsArg* leaf ;
624  while((leaf=(RooAbsArg*)liter->Next())) {
625  if (leaf->operMode()==ADirty && leaf->isValueServer(*this)) {
627  break ;
628  }
629  if (leaf->getAttribute("projectedDependent")) {
631  break ;
632  }
633  }
634  delete liter ;
635  }
636  }
637  delete siter ;
638 }
639 
640 
641 
642 ////////////////////////////////////////////////////////////////////////////////
643 /// Utility function that returns true if 'object server' is a server
644 /// to exactly one of the RooAbsArgs in 'exclLVBranches'
645 
646 Bool_t RooRealIntegral::servesExclusively(const RooAbsArg* server,const RooArgSet& exclLVBranches, const RooArgSet& allBranches) const
647 {
648  // Determine if given server serves exclusively exactly one of the given nodes in exclLVBranches
649 
650  // Special case, no LV servers available
651  if (exclLVBranches.getSize()==0) return kFALSE ;
652 
653  // If server has no clients and is not an LValue itself, return false
654  if (server->_clientList.GetSize()==0 && exclLVBranches.find(server->GetName())) {
655  return kFALSE ;
656  }
657 
658  // WVE must check for value relations only here!!!!
659 
660 
661 // cout << "servesExclusively: does " << server->GetName() << " serve only one of " << exclLVBranches << endl ;
662 
663  // Loop over all clients
664  Int_t numLVServ(0) ;
665  RooAbsArg* client ;
666  TIterator* cIter = server->valueClientIterator() ;
667  while((client=(RooAbsArg*)cIter->Next())) {
668 // cout << "now checking value client " << client->GetName() << " of server " << server->GetName() << endl ;
669  // If client is not an LValue, recurse
670  if (!(exclLVBranches.find(client->GetName())==client)) {
671 // cout << " client " << client->GetName() << "is not an lvalue" << endl ;
672  if (allBranches.find(client->GetName())==client) {
673 // cout << " ... recursing call" << endl ;
674  if (!servesExclusively(client,exclLVBranches,allBranches)) {
675  // Client is a non-LValue that doesn't have an exclusive LValue server
676  delete cIter ;
677 // cout << "client " << client->GetName() << " is a non-lvalue that doesn't have an exclusive lvalue server" << endl ;
678  return kFALSE ;
679  }
680  }
681  } else {
682  // Client is an LValue
683 // cout << "client " << client->GetName() << " of server " << server->GetName() << " is an LValue " << endl ;
684  numLVServ++ ;
685  }
686  }
687 
688  delete cIter ;
689 // cout << "numLVserv = " << numLVServ << endl ;
690  return (numLVServ==1) ;
691 }
692 
693 
694 
695 
696 ////////////////////////////////////////////////////////////////////////////////
697 /// (Re)Initialize numerical integration engine if necessary. Return kTRUE if
698 /// successful, or otherwise kFALSE.
699 
701 {
702  // if we already have an engine, check if it still works for the present limits.
703  if(0 != _numIntEngine) {
705  // otherwise, cleanup the old engine
706  delete _numIntEngine ;
707  _numIntEngine= 0;
708  if(0 != _numIntegrand) {
709  delete _numIntegrand;
710  _numIntegrand= 0;
711  }
712  }
713 
714  // All done if there are no arguments to integrate numerically
715  if(0 == _intList.getSize()) return kTRUE;
716 
717  // Bind the appropriate analytic integral (specified by _mode) of our RooRealVar object to
718  // those of its arguments that will be integrated out numerically.
719  if(_mode != 0) {
721  }
722  else {
724  }
725  if(0 == _numIntegrand || !_numIntegrand->isValid()) {
726  coutE(Integration) << ClassName() << "::" << GetName() << ": failed to create valid integrand." << endl;
727  return kFALSE;
728  }
729 
730  // Create appropriate numeric integrator using factory
733 
734  if(0 == _numIntEngine || !_numIntEngine->isValid()) {
735  coutE(Integration) << ClassName() << "::" << GetName() << ": failed to create valid integrator." << endl;
736  return kFALSE;
737  }
738 
739  cxcoutI(NumIntegration) << "RooRealIntegral::init(" << GetName() << ") using numeric integrator "
740  << _numIntEngine->IsA()->GetName() << " to calculate Int" << _intList << endl ;
741 
742  if (_intList.getSize()>3) {
743  cxcoutI(NumIntegration) << "RooRealIntegral::init(" << GetName() << ") evaluation requires " << _intList.getSize() << "-D numeric integration step. Evaluation may be slow, sufficient numeric precision for fitting & minimization is not guaranteed" << endl ;
744  }
745 
747  return kTRUE;
748 }
749 
750 
751 
752 ////////////////////////////////////////////////////////////////////////////////
753 /// Copy constructor
754 
756  RooAbsReal(other,name),
757  _valid(other._valid),
758  _sumList("!sumList",this,other._sumList),
759  _intList("!intList",this,other._intList),
760  _anaList("!anaList",this,other._anaList),
761  _jacList("!jacList",this,other._jacList),
762  _facList("!facList","Variables independent of function",this,kFALSE,kTRUE),
763  _facListIter(_facList.createIterator()),
764  _jacListIter(_jacList.createIterator()),
765  _function("!func",this,other._function),
766  _iconfig(other._iconfig),
767  _sumCat("!sumCat",this,other._sumCat),
768  _sumCatIter(0),
769  _mode(other._mode),
770  _intOperMode(other._intOperMode),
772  _numIntEngine(0),
773  _numIntegrand(0),
774  _rangeName(other._rangeName),
775  _params(0),
777 {
779 
780  other._facListIter->Reset() ;
781  RooAbsArg* arg ;
782  while((arg=(RooAbsArg*)other._facListIter->Next())) {
783  RooAbsArg* argClone = (RooAbsArg*) arg->Clone() ;
784  _facListOwned.addOwned(*argClone) ;
785  _facList.add(*argClone) ;
786  addServer(*argClone,kFALSE,kTRUE) ;
787  }
788 
789  other._intList.snapshot(_saveInt) ;
790  other._sumList.snapshot(_saveSum) ;
791 
793 }
794 
795 
796 
797 ////////////////////////////////////////////////////////////////////////////////
798 
800  // Destructor
801 {
802  if (_numIntEngine) delete _numIntEngine ;
803  if (_numIntegrand) delete _numIntegrand ;
804  if (_funcNormSet) delete _funcNormSet ;
805  delete _facListIter ;
806  delete _jacListIter ;
807  if (_sumCatIter) delete _sumCatIter ;
808  if (_params) delete _params ;
809 
811 }
812 
813 
814 
815 
816 
817 ////////////////////////////////////////////////////////////////////////////////
818 
819 RooAbsReal* RooRealIntegral::createIntegral(const RooArgSet& iset, const RooArgSet* nset, const RooNumIntConfig* cfg, const char* rangeName) const
820 {
821  // Handle special case of no integration with default algorithm
822  if (iset.getSize()==0) {
823  return RooAbsReal::createIntegral(iset,nset,cfg,rangeName) ;
824  }
825 
826  // Special handling of integral of integral, return RooRealIntegral that represents integral over all dimensions in one pass
827  RooArgSet isetAll(iset) ;
828  isetAll.add(_sumList) ;
829  isetAll.add(_intList) ;
830  isetAll.add(_anaList) ;
831  isetAll.add(_facList) ;
832 
833  const RooArgSet* newNormSet(0) ;
834  RooArgSet* tmp(0) ;
835  if (nset && !_funcNormSet) {
836  newNormSet = nset ;
837  } else if (!nset && _funcNormSet) {
838  newNormSet = _funcNormSet ;
839  } else if (nset && _funcNormSet) {
840  tmp = new RooArgSet ;
841  tmp->add(*nset) ;
842  tmp->add(*_funcNormSet,kTRUE) ;
843  newNormSet = tmp ;
844  }
845  RooAbsReal* ret = _function.arg().createIntegral(isetAll,newNormSet,cfg,rangeName) ;
846 
847  if (tmp) {
848  delete tmp ;
849  }
850 
851  return ret ;
852 }
853 
854 
855 
856 ////////////////////////////////////////////////////////////////////////////////
857 /// Return value of object. If the cache is clean, return the
858 /// cached value, otherwise recalculate on the fly and refill
859 /// the cache
860 
862 {
863 // // fast-track clean-cache processing
864 // if (_operMode==AClean) {
865 // return _value ;
866 // }
867 
868  if (nset && nset!=_lastNSet) {
869  ((RooAbsReal*) this)->setProxyNormSet(nset) ;
870  _lastNSet = (RooArgSet*) nset ;
871  }
872 
874  _value = traceEval(nset) ;
875  }
876 
877  return _value ;
878 }
879 
880 
881 
882 
883 
884 ////////////////////////////////////////////////////////////////////////////////
885 /// Perform the integration and return the result
886 
888 {
889  Double_t retVal(0) ;
890  switch (_intOperMode) {
891 
892  case Hybrid:
893  {
894  // Cache numeric integrals in >1d expensive object cache
895  RooDouble* cacheVal(0) ;
898  }
899 
900  if (cacheVal) {
901  retVal = *cacheVal ;
902  // cout << "using cached value of integral" << GetName() << endl ;
903  } else {
904 
905 
906  // Find any function dependents that are AClean
907  // and switch them temporarily to ADirty
908  Bool_t origState = inhibitDirty() ;
910 
911  // try to initialize our numerical integration engine
912  if(!(_valid= initNumIntegrator())) {
913  coutE(Integration) << ClassName() << "::" << GetName()
914  << ":evaluate: cannot initialize numerical integrator" << endl;
915  return 0;
916  }
917 
918  // Save current integral dependent values
919  _saveInt = _intList ;
920  _saveSum = _sumList ;
921 
922  // Evaluate sum/integral
923  retVal = sum() ;
924 
925  // This must happen BEFORE restoring dependents, otherwise no dirty state propagation in restore step
926  setDirtyInhibit(origState) ;
927 
928  // Restore integral dependent values
931 
932  // Cache numeric integrals in >1d expensive object cache
934  RooDouble* val = new RooDouble(retVal) ;
936 // cout << "### caching value of integral" << GetName() << " in " << &expensiveObjectCache() << endl ;
937  }
938 
939  }
940  break ;
941  }
942  case Analytic:
943  {
945  cxcoutD(Tracing) << "RooRealIntegral::evaluate_analytic(" << GetName()
946  << ")func = " << _function.arg().IsA()->GetName() << "::" << _function.arg().GetName()
947  << " raw = " << retVal << " _funcNormSet = " << (_funcNormSet?*_funcNormSet:RooArgSet()) << endl ;
948 
949 
950  break ;
951  }
952 
953  case PassThrough:
954  {
955  //setDirtyInhibit(kTRUE) ;
956  retVal= _function.arg().getVal(_funcNormSet) ;
957  //setDirtyInhibit(kFALSE) ;
958  break ;
959  }
960  }
961 
962 
963  // Multiply answer with integration ranges of factorized variables
964  if (_facList.getSize()>0) {
965  RooAbsArg *arg ;
966  _facListIter->Reset() ;
967  while((arg=(RooAbsArg*)_facListIter->Next())) {
968  // Multiply by fit range for 'real' dependents
969  if (arg->IsA()->InheritsFrom(RooAbsRealLValue::Class())) {
970  RooAbsRealLValue* argLV = (RooAbsRealLValue*)arg ;
971  retVal *= (argLV->getMax() - argLV->getMin()) ;
972  }
973  // Multiply by number of states for category dependents
974  if (arg->IsA()->InheritsFrom(RooAbsCategoryLValue::Class())) {
976  retVal *= argLV->numTypes() ;
977  }
978  }
979  }
980 
981 
982  if (dologD(Tracing)) {
983  cxcoutD(Tracing) << "RooRealIntegral::evaluate(" << GetName() << ") anaInt = " << _anaList << " numInt = " << _intList << _sumList << " mode = " ;
984  switch(_mode) {
985  case Hybrid: ccoutD(Tracing) << "Hybrid" ; break ;
986  case Analytic: ccoutD(Tracing) << "Analytic" ; break ;
987  case PassThrough: ccoutD(Tracing) << "PassThrough" ; break ;
988  }
989 
990  ccxcoutD(Tracing) << "raw*fact = " << retVal << endl ;
991  }
992 
993 
994  return retVal ;
995 }
996 
997 
998 
999 ////////////////////////////////////////////////////////////////////////////////
1000 /// Return product of jacobian terms originating from analytical integration
1001 
1003 {
1004  if (_jacList.getSize()==0) {
1005  return 1 ;
1006  }
1007 
1008  Double_t jacProd(1) ;
1009  _jacListIter->Reset() ;
1010  RooAbsRealLValue* arg ;
1011  while ((arg=(RooAbsRealLValue*)_jacListIter->Next())) {
1012  jacProd *= arg->jacobian() ;
1013  }
1014 
1015  // Take fabs() here: if jacobian is negative, min and max are swapped and analytical integral
1016  // will be positive, so must multiply with positive jacobian.
1017  return fabs(jacProd) ;
1018 }
1019 
1020 
1021 
1022 ////////////////////////////////////////////////////////////////////////////////
1023 /// Perform summation of list of category dependents to be integrated
1024 
1026 {
1027  if (_sumList.getSize()!=0) {
1028 
1029  // Add integrals for all permutations of categories summed over
1030  Double_t total(0) ;
1031 
1032  _sumCatIter->Reset() ;
1033  RooCatType* type ;
1035  while((type=(RooCatType*)_sumCatIter->Next())) {
1036  sumCat->setIndex(type->getVal()) ;
1037  if (!_rangeName || sumCat->inRange(RooNameReg::str(_rangeName))) {
1038  total += integrate() / jacobianProduct() ;
1039  }
1040  }
1041 
1042  return total ;
1043 
1044  } else {
1045 
1046  // Simply return integral
1047  Double_t ret = integrate() / jacobianProduct() ;
1048  return ret ;
1049  }
1050 }
1051 
1052 
1053 
1054 
1055 ////////////////////////////////////////////////////////////////////////////////
1056 /// Perform hybrid numerical/analytical integration over all real-valued dependents
1057 
1059 {
1060  if (!_numIntEngine) {
1061  // Trivial case, fully analytical integration
1063  }
1064  else {
1065  return _numIntEngine->calculate() ;
1066  }
1067 }
1068 
1069 
1070 
1071 ////////////////////////////////////////////////////////////////////////////////
1072 /// Intercept server redirects and reconfigure internal object accordingly
1073 
1075  Bool_t /*mustReplaceAll*/, Bool_t /*nameChange*/, Bool_t /*isRecursive*/)
1076 {
1078 
1080 
1081  // Update contents value caches for _intList and _sumList
1082  _saveInt.removeAll() ;
1083  _saveSum.removeAll() ;
1086 
1087  // Delete parameters cache if we have one
1088  if (_params) {
1089  delete _params ;
1090  _params = 0 ;
1091  }
1092 
1093  return kFALSE ;
1094 }
1095 
1096 
1097 
1098 ////////////////////////////////////////////////////////////////////////////////
1099 
1101 {
1102  if (!_params) {
1103  _params = new RooArgSet("params") ;
1104 
1105  TIterator* siter = serverIterator() ;
1106  RooArgSet params ;
1107  RooAbsArg* server ;
1108  while((server = (RooAbsArg*)siter->Next())) {
1109  if (server->isValueServer(*this)) _params->add(*server) ;
1110  }
1111  delete siter ;
1112  }
1113 
1114  return *_params ;
1115 }
1116 
1117 
1118 
1119 ////////////////////////////////////////////////////////////////////////////////
1120 /// Dummy
1121 
1123 {
1124  if (_operMode==ADirty) {
1125 // cout << "RooRealIntegral::operModeHook(" << GetName() << " warning: mode set to ADirty" << endl ;
1126 // if (TString(GetName()).Contains("FULL")) {
1127 // cout << "blah" << endl ;
1128 // }
1129  }
1130 }
1131 
1132 
1133 
1134 ////////////////////////////////////////////////////////////////////////////////
1135 /// Check if current value is valid
1136 
1137 Bool_t RooRealIntegral::isValidReal(Double_t /*value*/, Bool_t /*printError*/) const
1138 {
1139  return kTRUE ;
1140 }
1141 
1142 
1143 
1144 ////////////////////////////////////////////////////////////////////////////////
1145 /// Customized printing of arguments of a RooRealIntegral to more intuitively reflect the contents of the
1146 /// integration operation
1147 
1148 void RooRealIntegral::printMetaArgs(ostream& os) const
1149 {
1150 
1151  if (intVars().getSize()!=0) {
1152  os << "Int " ;
1153  }
1154  os << _function.arg().GetName() ;
1155  if (_funcNormSet) {
1156  os << "_Norm" ;
1157  os << *_funcNormSet ;
1158  os << " " ;
1159  }
1160 
1161  // List internally integrated observables and factorizing observables as analytically integrated
1162  RooArgSet tmp(_anaList) ;
1163  tmp.add(_facList) ;
1164  if (tmp.getSize()>0) {
1165  os << "d[Ana]" ;
1166  os << tmp ;
1167  os << " " ;
1168  }
1169 
1170  // List numerically integrated and summed observables as numerically integrated
1171  RooArgSet tmp2(_intList) ;
1172  tmp2.add(_sumList) ;
1173  if (tmp2.getSize()>0) {
1174  os << " d[Num]" ;
1175  os << tmp2 ;
1176  os << " " ;
1177  }
1178 }
1179 
1180 
1181 
1182 ////////////////////////////////////////////////////////////////////////////////
1183 /// Print the state of this object to the specified output stream.
1184 
1185 void RooRealIntegral::printMultiline(ostream& os, Int_t contents, Bool_t verbose, TString indent) const
1186 {
1187  RooAbsReal::printMultiline(os,contents,verbose,indent) ;
1188  os << indent << "--- RooRealIntegral ---" << endl;
1189  os << indent << " Integrates ";
1191  TString deeper(indent);
1192  deeper.Append(" ");
1193  os << indent << " operating mode is "
1194  << (_intOperMode==Hybrid?"Hybrid":(_intOperMode==Analytic?"Analytic":"PassThrough")) << endl ;
1195  os << indent << " Summed discrete args are " << _sumList << endl ;
1196  os << indent << " Numerically integrated args are " << _intList << endl;
1197  os << indent << " Analytically integrated args using mode " << _mode << " are " << _anaList << endl ;
1198  os << indent << " Arguments included in Jacobian are " << _jacList << endl ;
1199  os << indent << " Factorized arguments are " << _facList << endl ;
1200  os << indent << " Function normalization set " ;
1201  if (_funcNormSet)
1202  _funcNormSet->Print("1") ;
1203  else
1204  os << "<none>" ;
1205 
1206  os << endl ;
1207 }
1208 
1209 
1210 
1211 ////////////////////////////////////////////////////////////////////////////////
1212 /// Global switch to cache all integral values that integrate at least ndim dimensions numerically
1213 
1215  _cacheAllNDim = ndim ;
1216 }
1217 
1218 
1219 ////////////////////////////////////////////////////////////////////////////////
1220 /// Return minimum dimensions of numeric integration for which values are cached.
1221 
1223 {
1224  return _cacheAllNDim ;
1225 }
1226 
1227 
virtual Double_t getMin(const char *name=0) const
RooSetProxy _anaList
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:51
RooSetProxy _sumList
TIterator * createIterator(Bool_t dir=kIterForward) const
Bool_t initNumIntegrator() const
(Re)Initialize numerical integration engine if necessary.
#define coutE(a)
Definition: RooMsgService.h:35
virtual void printStream(std::ostream &os, Int_t contents, StyleOption style, TString indent="") const
Print description of object on ostream, printing contents set by contents integer, which is interpreted as an OR of &#39;enum ContentsOptions&#39; values and in the style given by &#39;enum StyleOption&#39;.
RooArgSet intVars() const
RooListProxy _sumCat
virtual ~RooRealIntegral()
RooNumIntConfig holds the configuration parameters of the various numeric integrators used by RooReal...
virtual Double_t getMax(const char *name=0) const
RooRealProxy _function
Iterator over lvalue observables with Jacobian.
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:86
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:744
TIterator * _facListIter
#define cxcoutI(a)
Definition: RooMsgService.h:84
virtual void Reset()=0
virtual Bool_t isValidReal(Double_t value, Bool_t printError=kFALSE) const
Check if current value is valid.
static RooNumIntFactory & instance()
Static method returning reference to singleton instance of factory.
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Definition: RooAbsArg.h:194
void leafNodeServerList(RooAbsCollection *list, const RooAbsArg *arg=0, Bool_t recurseNonDerived=kFALSE) const
Fill supplied list with all leaf nodes of the arg tree, starting with ourself as top node...
Definition: RooAbsArg.cxx:493
virtual Bool_t isParameterized() const
Definition: RooAbsBinning.h:79
virtual void removeAll()
Remove all argument inset using remove(const RooAbsArg&).
Double_t getVal(const RooArgSet *set=0) const
Definition: RooAbsReal.h:64
#define cxcoutD(a)
Definition: RooMsgService.h:80
static const char * str(const TNamed *ptr)
Return C++ string corresponding to given TNamed pointer.
Definition: RooNameReg.cxx:135
Bool_t isValueOrShapeDirtyAndClear() const
Definition: RooAbsArg.h:370
Double_t traceEval(const RooArgSet *set) const
Calculate current value of object, with error tracing wrapper.
Definition: RooAbsReal.cxx:288
Bool_t _restartNumIntEngine
virtual void printMultiline(std::ostream &os, Int_t contents, Bool_t verbose=kFALSE, TString indent="") const
Structure printing.
Definition: RooAbsReal.cxx:424
virtual Double_t jacobianProduct() const
Return product of jacobian terms originating from analytical integration.
virtual Bool_t isLValue() const
Definition: RooAbsArg.h:170
int Int_t
Definition: RtypesCore.h:41
void printMultiline(std::ostream &os, Int_t contents, Bool_t verbose=kFALSE, TString indent="") const
Print the state of this object to the specified output stream.
bool Bool_t
Definition: RtypesCore.h:59
Double_t calculate(const Double_t *yvec=0)
Calculate integral value with given array of parameter values.
TIterator * valueClientIterator() const
Definition: RooAbsArg.h:104
const Bool_t kFALSE
Definition: Rtypes.h:92
Double_t evaluate() const
Perform the integration and return the result.
virtual TObject * Clone(const char *newname=0) const
Make a clone of an object using the Streamer facility.
Definition: RooAbsArg.h:75
const RooArgSet & parameters() const
TIterator * _jacListIter
Iterator over factorizing observables.
Bool_t inhibitDirty() const
Delete watch flag.
Definition: RooAbsArg.cxx:84
RooArgSet * _lastNSet
Definition: RooAbsReal.h:481
RooArgSet * _params
STL namespace.
const char * Class
Definition: TXMLSetup.cxx:64
Bool_t isValid() const
Bool_t isValid() const
Definition: RooAbsFunc.h:33
virtual Double_t jacobian() const
RooDouble is a minimal implementation of a TObject holding a Double_t value.
Definition: RooDouble.h:22
const TObject * retrieveObject(const char *name, TClass *tclass, const RooArgSet &params)
Retrieve object from cache that was registered under given name with given parameters, if current parameter values match those that were stored in the registry for this object.
Int_t GetSize() const
Definition: RooLinkedList.h:60
Iterator abstract base class.
Definition: TIterator.h:32
RooSetProxy _intList
friend class RooRealBinding
Definition: RooAbsReal.h:385
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition: TObject.cxx:188
#define TRACE_CREATE
Definition: RooTrace.h:23
RooRealIntegral performs hybrid numerical/analytical integrals of RooAbsReal objects The class perfor...
Int_t numTypes(const char *=0) const
RooAbsIntegrator * createIntegrator(RooAbsFunc &func, const RooNumIntConfig &config, Int_t ndim=0, Bool_t isBinned=kFALSE)
Construct a numeric integrator instance that operates on function &#39;func&#39; and is configured with &#39;conf...
RooAbsIntegrator * _numIntEngine
do not persist
virtual void removeAll()
Remove all arguments from our set, deleting them if we own them.
RooNumIntConfig * _iconfig
do not persist
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:33
virtual Bool_t checkLimits() const
virtual void Print(Option_t *options=0) const
Print TNamed name and title.
Definition: RooAbsArg.h:227
if on multiple lines(like in C++). **The " * configuration fragment. * * The "import myobject continue
Parses the configuration file.
Definition: HLFactory.cxx:368
#define oocxcoutD(o, a)
Definition: RooMsgService.h:82
#define ccxcoutD(a)
Definition: RooMsgService.h:81
RooSetProxy _jacList
virtual RooAbsReal * highBoundFunc() const
Definition: RooAbsBinning.h:87
RooCatType is an auxilary class for RooAbsCategory and defines a a single category state...
Definition: RooCatType.h:23
friend class RooArgSet
Definition: RooAbsArg.h:469
virtual RooAbsReal * lowBoundFunc() const
Definition: RooAbsBinning.h:83
virtual void Print(Option_t *options=0) const
This method must be overridden when a class wants to print itself.
RooAbsCategoryLValue is the common abstract base class for objects that represent a discrete value th...
Bool_t overlaps(const RooAbsArg &testArg, Bool_t valueOnly=kFALSE) const
Test if any of the nodes of tree are shared with that of the given tree.
Definition: RooAbsArg.cxx:803
virtual Double_t sum() const
Perform summation of list of category dependents to be integrated.
static Int_t getCacheAllNumeric()
Return minimum dimensions of numeric integration for which values are cached.
virtual void addClone(const RooAbsCollection &col, Bool_t silent=kFALSE)
Add a collection of arguments to this collection by calling addOwned() for each element in the source...
Definition: RooArgSet.h:94
virtual Bool_t isFundamental() const
Definition: RooAbsArg.h:157
Bool_t getAttribute(const Text_t *name) const
Check if a named attribute is set. By default, all attributes are unset.
Definition: RooAbsArg.cxx:289
RooArgSet _saveSum
do not persist
Int_t getSize() const
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:503
VecExpr< UnaryOp< Fabs< T >, VecExpr< A, T, D >, T >, T, D > fabs(const VecExpr< A, T, D > &rhs)
RooAbsCollection * snapshot(Bool_t deepCopy=kTRUE) const
Take a snap shot of current collection contents: An owning collection is returned containing clones o...
RooExpensiveObjectCache & expensiveObjectCache() const
Definition: RooAbsArg.cxx:2342
TIterator * typeIterator() const
Return iterator over all defined states.
RooAbsArg * first() const
RooSetProxy _facList
RooArgSet _facListOwned
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:488
virtual void operModeHook()
Cache all integrals with given numeric dimension.
virtual Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &analVars, const RooArgSet *normSet, const char *rangeName=0) const
Variant of getAnalyticalIntegral that is also passed the normalization set that should be applied to ...
Definition: RooAbsReal.cxx:320
bool verbose
char * Form(const char *fmt,...)
RooRefCountList _clientList
Definition: RooAbsArg.h:474
virtual Bool_t addOwned(const RooAbsCollection &col, Bool_t silent=kFALSE)
Add a collection of arguments to this collection by calling addOwned() for each element in the source...
Definition: RooArgSet.h:90
virtual const RooAbsBinning & getBinning(const char *name=0, Bool_t verbose=kTRUE, Bool_t createOnTheFly=kFALSE) const =0
OperMode operMode() const
Definition: RooAbsArg.h:399
OperMode _operMode
Definition: RooAbsArg.h:570
TIterator * _sumCatIter
do not persist
virtual Bool_t setIndex(Int_t index, Bool_t printError=kTRUE)
Set the value of the super category by specifying the state index code by setting the states of the c...
static void indent(ostringstream &buf, int indent_level)
IntOperMode _intOperMode
TIterator * serverIterator() const
Definition: RooAbsArg.h:112
static void setDirtyInhibit(Bool_t flag)
Control global dirty inhibit mode.
Definition: RooAbsArg.cxx:239
virtual Bool_t forceAnalyticalInt(const RooAbsArg &) const
Definition: RooAbsReal.h:104
static Int_t _cacheAllNDim
static unsigned int total
Lightweight RooAbsFunc interface adaptor that binds an analytic integral of a RooAbsReal object (spec...
virtual void setExpensiveObjectCache(RooExpensiveObjectCache &cache)
Definition: RooAbsArg.h:497
void autoSelectDirtyMode()
Set appropriate cache operation mode for integral depending on cache operation mode of server objects...
virtual Bool_t isBinnedDistribution(const RooArgSet &) const
Definition: RooAbsReal.h:277
virtual Bool_t isJacobianOK(const RooArgSet &depList) const
#define ClassImp(name)
Definition: Rtypes.h:279
RooAbsArg * find(const char *name) const
Find object with given name in list.
Bool_t servesExclusively(const RooAbsArg *server, const RooArgSet &exclLVBranches, const RooArgSet &allBranches) const
Utility function that returns true if &#39;object server&#39; is a server to exactly one of the RooAbsArgs in...
double Double_t
Definition: RtypesCore.h:55
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:53
virtual Double_t integrate() const
Perform hybrid numerical/analytical integration over all real-valued dependents.
int type
Definition: TGX11.cxx:120
#define TRACE_DESTROY
Definition: RooTrace.h:24
virtual Bool_t addOwned(RooAbsArg &var, Bool_t silent=kFALSE)
Reimplementation of standard RooArgList::addOwned()
Bool_t _cacheNum
cache for set of parameters
const char * GetName() const
Returns name of object.
Double_t _value
Definition: RooAbsReal.h:389
Bool_t dependsOnValue(const RooAbsCollection &serverList, const RooAbsArg *ignoreArg=0) const
Definition: RooAbsArg.h:88
#define dologD(a)
Definition: RooMsgService.h:64
RooAbsCollection is an abstract container object that can hold multiple RooAbsArg objects...
static void setCacheAllNumeric(Int_t ndim)
Global switch to cache all integral values that integrate at least ndim dimensions numerically...
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE)
Remove the specified argument from our list.
Bool_t isShapeServer(const RooAbsArg &arg) const
Definition: RooAbsArg.h:142
virtual Bool_t redirectServersHook(const RooAbsCollection &newServerList, Bool_t mustReplaceAll, Bool_t nameChange, Bool_t isRecursive)
Intercept server redirects and reconfigure internal object accordingly.
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
Bool_t isValueServer(const RooAbsArg &arg) const
Definition: RooAbsArg.h:134
virtual Bool_t isDerived() const
Definition: RooAbsArg.h:81
virtual RooAbsReal * createIntegral(const RooArgSet &iset, const RooArgSet *nset=0, const RooNumIntConfig *cfg=0, const char *rangeName=0) const
Create an object that represents the integral of the function over one or more observables listed in ...
Int_t getVal() const
Definition: RooCatType.h:80
const RooNumIntConfig * getIntegratorConfig() const
Return the numeric integration configuration used for this object.
#define oocxcoutI(o, a)
Definition: RooMsgService.h:86
RooNameReg is a registry for &#39;const char*&#39; name.
Definition: RooNameReg.h:23
virtual TObject * Next()=0
RooArgSet * _funcNormSet
RooAbsFunc * _numIntegrand
do not persist
void setOperMode(OperMode mode, Bool_t recurseADirty=kTRUE)
Change cache operation mode to given mode.
Definition: RooAbsArg.cxx:1754
RooSuperCategory consolidates several RooAbsCategoryLValue objects into a single category.
Bool_t contains(const RooAbsArg &var) const
RooAbsArg is the common abstract base class for objects that represent a value (of arbitrary type) an...
Definition: RooAbsArg.h:66
const RooAbsReal & arg() const
Definition: RooRealProxy.h:43
void printMetaArgs(std::ostream &os) const
Customized printing of arguments of a RooRealIntegral to more intuitively reflect the contents of the...
virtual const char * GetName() const
Returns name of object.
Definition: TObject.cxx:416
void addServer(RooAbsArg &server, Bool_t valueProp=kTRUE, Bool_t shapeProp=kFALSE)
Register another RooAbsArg as a server to us, ie, declare that we depend on it.
Definition: RooAbsArg.cxx:362
TNamed * _rangeName
do not persist
const Bool_t kTRUE
Definition: Rtypes.h:91
#define ccoutD(a)
Definition: RooMsgService.h:38
virtual Bool_t inRange(const char *rangeName) const
Return true of all of the input category states are in the given range.
virtual Double_t getValV(const RooArgSet *set=0) const
Return value of object.
Bool_t registerObject(const char *ownerName, const char *objectName, TObject &cacheObject, TIterator *paramIter)
Register object associated with given name and given associated parameters with given values in cache...
char name[80]
Definition: TGX11.cxx:109
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Overloaded RooArgSet::add() method inserts &#39;var&#39; into set and registers &#39;var&#39; as server to owner with...
virtual Double_t analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=0) const
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral. ...
Definition: RooAbsReal.cxx:349