Logo ROOT   6.10/09
Reference Guide
RooSimultaneous.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 RooSimultaneous.cxx
19 \class RooSimultaneous
20 \ingroup Roofitcore
21 
22 RooSimultaneous facilitates simultaneous fitting of multiple PDFs
23 to subsets of a given dataset.
24 The class takes an index category, which is interpreted as
25 the data subset indicator, and a list of PDFs, each associated
26 with a state of the index category. RooSimultaneous always returns
27 the value of the PDF that is associated with the current value
28 of the index category
29 Extended likelihood fitting is supported if all components support
30 extended likelihood mode. The expected number of events by a RooSimultaneous
31 is that of the component p.d.f. selected by the index category
32 **/
33 
34 #include "RooFit.h"
35 #include "Riostream.h"
36 
37 #include "TObjString.h"
38 #include "RooSimultaneous.h"
39 #include "RooAbsCategoryLValue.h"
40 #include "RooPlot.h"
41 #include "RooCurve.h"
42 #include "RooRealVar.h"
43 #include "RooAddPdf.h"
44 #include "RooAbsData.h"
45 #include "Roo1DTable.h"
46 #include "RooSimGenContext.h"
47 #include "RooSimSplitGenContext.h"
48 #include "RooDataSet.h"
49 #include "RooCmdConfig.h"
50 #include "RooNameReg.h"
51 #include "RooGlobalFunc.h"
52 #include "RooNameReg.h"
53 #include "RooMsgService.h"
54 #include "RooCategory.h"
55 #include "RooSuperCategory.h"
56 #include "RooDataHist.h"
57 #include "RooRandom.h"
58 #include "RooArgSet.h"
59 
60 using namespace std ;
61 
63 ;
64 
65 
66 
67 
68 ////////////////////////////////////////////////////////////////////////////////
69 /// Constructor with index category. PDFs associated with indexCat
70 /// states can be added after construction with the addPdf() function.
71 ///
72 /// RooSimultaneous can function without having a PDF associated
73 /// with every single state. The normalization in such cases is taken
74 /// from the number of registered PDFs, but getVal() will assert if
75 /// when called for an unregistered index state.
76 
77 RooSimultaneous::RooSimultaneous(const char *name, const char *title,
78  RooAbsCategoryLValue& inIndexCat) :
79  RooAbsPdf(name,title),
80  _plotCoefNormSet("!plotCoefNormSet","plotCoefNormSet",this,kFALSE,kFALSE),
81  _plotCoefNormRange(0),
82  _partIntMgr(this,10),
83  _indexCat("indexCat","Index category",this,inIndexCat),
84  _numPdf(0)
85 {
86 }
87 
88 
89 
90 ////////////////////////////////////////////////////////////////////////////////
91 /// Constructor from index category and full list of PDFs.
92 /// In this constructor form, a PDF must be supplied for each indexCat state
93 /// to avoid ambiguities. The PDFS are associated in order with the state of the
94 /// index category as listed by the index categories type iterator.
95 ///
96 /// PDFs may not overlap (i.e. share any variables) with the index category (function)
97 
98 RooSimultaneous::RooSimultaneous(const char *name, const char *title,
99  const RooArgList& inPdfList, RooAbsCategoryLValue& inIndexCat) :
100  RooAbsPdf(name,title),
101  _plotCoefNormSet("!plotCoefNormSet","plotCoefNormSet",this,kFALSE,kFALSE),
103  _partIntMgr(this,10),
104  _indexCat("indexCat","Index category",this,inIndexCat),
105  _numPdf(0)
106 {
107  if (inPdfList.getSize() != inIndexCat.numTypes()) {
108  coutE(InputArguments) << "RooSimultaneous::ctor(" << GetName()
109  << " ERROR: Number PDF list entries must match number of index category states, no PDFs added" << endl ;
110  return ;
111  }
112 
113  map<string,RooAbsPdf*> pdfMap ;
114  // Iterator over PDFs and index cat states and add each pair
115  TIterator* pIter = inPdfList.createIterator() ;
116  TIterator* cIter = inIndexCat.typeIterator() ;
117  RooAbsPdf* pdf ;
118  RooCatType* type(0) ;
119  while ((pdf=(RooAbsPdf*)pIter->Next())) {
120  type = (RooCatType*) cIter->Next() ;
121  pdfMap[string(type->GetName())] = pdf ;
122  }
123  delete pIter ;
124  delete cIter ;
125 
126  initialize(inIndexCat,pdfMap) ;
127 }
128 
129 
130 ////////////////////////////////////////////////////////////////////////////////
131 
132 RooSimultaneous::RooSimultaneous(const char *name, const char *title,
133  map<string,RooAbsPdf*> pdfMap, RooAbsCategoryLValue& inIndexCat) :
134  RooAbsPdf(name,title),
135  _plotCoefNormSet("!plotCoefNormSet","plotCoefNormSet",this,kFALSE,kFALSE),
137  _partIntMgr(this,10),
138  _indexCat("indexCat","Index category",this,inIndexCat),
139  _numPdf(0)
140 {
141  initialize(inIndexCat,pdfMap) ;
142 }
143 
144 
145 
146 
147 // This class cannot be locally defined in initialize as it cannot be
148 // used as a template argument in that case
150  struct CompInfo {
151  RooAbsPdf* pdf ;
152  RooSimultaneous* simPdf ;
153  const RooAbsCategoryLValue* subIndex ;
154  RooArgSet* subIndexComps ;
155  } ;
156 }
157 
158 void RooSimultaneous::initialize(RooAbsCategoryLValue& inIndexCat, std::map<std::string,RooAbsPdf*> pdfMap)
159 {
160  // First see if there are any RooSimultaneous input components
161  Bool_t simComps(kFALSE) ;
162  for (map<string,RooAbsPdf*>::iterator iter=pdfMap.begin() ; iter!=pdfMap.end() ; iter++) {
163  if (dynamic_cast<RooSimultaneous*>(iter->second)) {
164  simComps = kTRUE ;
165  break ;
166  }
167  }
168 
169  // If there are no simultaneous component p.d.f. do simple processing through addPdf()
170  if (!simComps) {
171  for (map<string,RooAbsPdf*>::iterator iter=pdfMap.begin() ; iter!=pdfMap.end() ; iter++) {
172  addPdf(*iter->second,iter->first.c_str()) ;
173  }
174  return ;
175  }
176 
177  // Issue info message that we are about to do some rearraning
178  coutI(InputArguments) << "RooSimultaneous::initialize(" << GetName() << ") INFO: one or more input component of simultaneous p.d.f.s are"
179  << " simultaneous p.d.f.s themselves, rewriting composite expressions as one-level simultaneous p.d.f. in terms of"
180  << " final constituents and extended index category" << endl ;
181 
182 
183  RooArgSet allAuxCats ;
184  map<string,RooSimultaneousAux::CompInfo> compMap ;
185  for (map<string,RooAbsPdf*>::iterator iter=pdfMap.begin() ; iter!=pdfMap.end() ; iter++) {
186  RooSimultaneousAux::CompInfo ci ;
187  ci.pdf = iter->second ;
188  RooSimultaneous* simComp = dynamic_cast<RooSimultaneous*>(iter->second) ;
189  if (simComp) {
190  ci.simPdf = simComp ;
191  ci.subIndex = &simComp->indexCat() ;
192  ci.subIndexComps = simComp->indexCat().isFundamental() ? new RooArgSet(simComp->indexCat()) : simComp->indexCat().getVariables() ;
193  allAuxCats.add(*(ci.subIndexComps),kTRUE) ;
194  } else {
195  ci.simPdf = 0 ;
196  ci.subIndex = 0 ;
197  ci.subIndexComps = 0 ;
198  }
199  compMap[iter->first] = ci ;
200  }
201 
202  // Construct the 'superIndex' from the nominal index category and all auxiliary components
203  RooArgSet allCats(inIndexCat) ;
204  allCats.add(allAuxCats) ;
205  string siname = Form("%s_index",GetName()) ;
206  RooSuperCategory* superIndex = new RooSuperCategory(siname.c_str(),siname.c_str(),allCats) ;
207 
208  // Now process each of original pdf/state map entries
209  for (map<string,RooSimultaneousAux::CompInfo>::iterator citer = compMap.begin() ; citer != compMap.end() ; citer++) {
210 
211  RooArgSet repliCats(allAuxCats) ;
212  if (citer->second.subIndexComps) {
213  repliCats.remove(*citer->second.subIndexComps) ;
214  delete citer->second.subIndexComps ;
215  }
216  inIndexCat.setLabel(citer->first.c_str()) ;
217 
218 
219  if (!citer->second.simPdf) {
220 
221  // Entry is a plain p.d.f. assign it to every state permutation of the repliCats set
222  RooSuperCategory repliSuperCat("tmp","tmp",repliCats) ;
223 
224  // Iterator over all states of repliSuperCat
225  TIterator* titer = repliSuperCat.typeIterator() ;
226  RooCatType* type ;
227  while ((type=(RooCatType*)titer->Next())) {
228  // Set value
229  repliSuperCat.setLabel(type->GetName()) ;
230  // Retrieve corresponding label of superIndex
231  string superLabel = superIndex->getLabel() ;
232  addPdf(*citer->second.pdf,superLabel.c_str()) ;
233  cxcoutD(InputArguments) << "RooSimultaneous::initialize(" << GetName()
234  << ") assigning pdf " << citer->second.pdf->GetName() << " to super label " << superLabel << endl ;
235  }
236  } else {
237 
238  // Entry is a simultaneous p.d.f
239 
240  if (repliCats.getSize()==0) {
241 
242  // Case 1 -- No replication of components of RooSim component are required
243 
244  TIterator* titer = citer->second.subIndex->typeIterator() ;
245  RooCatType* type ;
246  while ((type=(RooCatType*)titer->Next())) {
247  const_cast<RooAbsCategoryLValue*>(citer->second.subIndex)->setLabel(type->GetName()) ;
248  string superLabel = superIndex->getLabel() ;
249  RooAbsPdf* compPdf = citer->second.simPdf->getPdf(type->GetName()) ;
250  if (compPdf) {
251  addPdf(*compPdf,superLabel.c_str()) ;
252  cxcoutD(InputArguments) << "RooSimultaneous::initialize(" << GetName()
253  << ") assigning pdf " << compPdf->GetName() << "(member of " << citer->second.pdf->GetName()
254  << ") to super label " << superLabel << endl ;
255  } else {
256  coutW(InputArguments) << "RooSimultaneous::initialize(" << GetName() << ") WARNING: No p.d.f. associated with label "
257  << type->GetName() << " for component RooSimultaneous p.d.f " << citer->second.pdf->GetName()
258  << "which is associated with master index label " << citer->first << endl ;
259  }
260  }
261  delete titer ;
262 
263  } else {
264 
265  // Case 2 -- Replication of components of RooSim component are required
266 
267  // Make replication supercat
268  RooSuperCategory repliSuperCat("tmp","tmp",repliCats) ;
269  TIterator* triter = repliSuperCat.typeIterator() ;
270 
271  TIterator* tsiter = citer->second.subIndex->typeIterator() ;
272  RooCatType* stype, *rtype ;
273  while ((stype=(RooCatType*)tsiter->Next())) {
274  const_cast<RooAbsCategoryLValue*>(citer->second.subIndex)->setLabel(stype->GetName()) ;
275  triter->Reset() ;
276  while ((rtype=(RooCatType*)triter->Next())) {
277  repliSuperCat.setLabel(rtype->GetName()) ;
278  string superLabel = superIndex->getLabel() ;
279  RooAbsPdf* compPdf = citer->second.simPdf->getPdf(stype->GetName()) ;
280  if (compPdf) {
281  addPdf(*compPdf,superLabel.c_str()) ;
282  cxcoutD(InputArguments) << "RooSimultaneous::initialize(" << GetName()
283  << ") assigning pdf " << compPdf->GetName() << "(member of " << citer->second.pdf->GetName()
284  << ") to super label " << superLabel << endl ;
285  } else {
286  coutW(InputArguments) << "RooSimultaneous::initialize(" << GetName() << ") WARNING: No p.d.f. associated with label "
287  << stype->GetName() << " for component RooSimultaneous p.d.f " << citer->second.pdf->GetName()
288  << "which is associated with master index label " << citer->first << endl ;
289  }
290  }
291  }
292 
293  delete tsiter ;
294  delete triter ;
295 
296  }
297  }
298  }
299 
300  // Change original master index to super index and take ownership of it
301  _indexCat.setArg(*superIndex) ;
302  addOwnedComponents(*superIndex) ;
303 
304 }
305 
306 
307 
308 ////////////////////////////////////////////////////////////////////////////////
309 /// Copy constructor
310 
312  RooAbsPdf(other,name),
313  _plotCoefNormSet("!plotCoefNormSet",this,other._plotCoefNormSet),
315  _partIntMgr(other._partIntMgr,this),
316  _indexCat("indexCat",this,other._indexCat),
317  _numPdf(other._numPdf)
318 {
319  // Copy proxy list
320  TIterator* pIter = other._pdfProxyList.MakeIterator() ;
321  RooRealProxy* proxy ;
322  while ((proxy=(RooRealProxy*)pIter->Next())) {
323  _pdfProxyList.Add(new RooRealProxy(proxy->GetName(),this,*proxy)) ;
324  }
325  delete pIter ;
326 }
327 
328 
329 
330 ////////////////////////////////////////////////////////////////////////////////
331 /// Destructor
332 
334 {
336 }
337 
338 
339 
340 ////////////////////////////////////////////////////////////////////////////////
341 /// Return the p.d.f associated with the given index category name
342 
343 RooAbsPdf* RooSimultaneous::getPdf(const char* catName) const
344 {
345  RooRealProxy* proxy = (RooRealProxy*) _pdfProxyList.FindObject(catName) ;
346  return proxy ? ((RooAbsPdf*)proxy->absArg()) : 0 ;
347 }
348 
349 
350 
351 ////////////////////////////////////////////////////////////////////////////////
352 /// Associate given PDF with index category state label 'catLabel'.
353 /// The names state must be already defined in the index category
354 ///
355 /// RooSimultaneous can function without having a PDF associated
356 /// with every single state. The normalization in such cases is taken
357 /// from the number of registered PDFs, but getVal() will assert if
358 /// when called for an unregistered index state.
359 ///
360 /// PDFs may not overlap (i.e. share any variables) with the index category (function)
361 
362 Bool_t RooSimultaneous::addPdf(const RooAbsPdf& pdf, const char* catLabel)
363 {
364  // PDFs cannot overlap with the index category
365  if (pdf.dependsOn(_indexCat.arg())) {
366  coutE(InputArguments) << "RooSimultaneous::addPdf(" << GetName() << "): ERROR, PDF " << pdf.GetName()
367  << " overlaps with index category " << _indexCat.arg().GetName() << endl ;
368  return kTRUE ;
369  }
370 
371  // Each index state can only have one PDF associated with it
372  if (_pdfProxyList.FindObject(catLabel)) {
373  coutE(InputArguments) << "RooSimultaneous::addPdf(" << GetName() << "): ERROR, index state "
374  << catLabel << " has already an associated PDF" << endl ;
375  return kTRUE ;
376  }
377 
378  const RooSimultaneous* simPdf = dynamic_cast<const RooSimultaneous*>(&pdf) ;
379  if (simPdf) {
380 
381  coutE(InputArguments) << "RooSimultaneous::addPdf(" << GetName()
382  << ") ERROR: you cannot add a RooSimultaneous component to a RooSimultaneous using addPdf()."
383  << " Use the constructor with RooArgList if input p.d.f.s or the map<string,RooAbsPdf&> instead." << endl ;
384  return kTRUE ;
385 
386  } else {
387 
388  // Create a proxy named after the associated index state
389  TObject* proxy = new RooRealProxy(catLabel,catLabel,this,(RooAbsPdf&)pdf) ;
390  _pdfProxyList.Add(proxy) ;
391  _numPdf += 1 ;
392  }
393 
394  return kFALSE ;
395 }
396 
397 
398 
399 
400 
401 ////////////////////////////////////////////////////////////////////////////////
402 /// WVE NEEDS FIX
403 
405 {
406  Bool_t allCanExtend(kTRUE) ;
407  Bool_t anyMustExtend(kFALSE) ;
408 
409  for (Int_t i=0 ; i<_numPdf ; i++) {
411  if (proxy) {
412 // cout << " now processing pdf " << pdf->GetName() << endl ;
413  RooAbsPdf* pdf = (RooAbsPdf*) proxy->absArg() ;
414  if (!pdf->canBeExtended()) {
415 // cout << "RooSim::extendedMode(" << GetName() << ") component " << pdf->GetName() << " cannot be extended" << endl ;
416  allCanExtend=kFALSE ;
417  }
418  if (pdf->mustBeExtended()) {
419  anyMustExtend=kTRUE;
420  }
421  }
422  }
423  if (anyMustExtend) {
424 // cout << "RooSim::extendedMode(" << GetName() << ") returning MustBeExtended" << endl ;
425  return MustBeExtended ;
426  }
427  if (allCanExtend) {
428 // cout << "RooSim::extendedMode(" << GetName() << ") returning CanBeExtended" << endl ;
429  return CanBeExtended ;
430  }
431 // cout << "RooSim::extendedMode(" << GetName() << ") returning CanNotBeExtended" << endl ;
432  return CanNotBeExtended ;
433 }
434 
435 
436 
437 
438 ////////////////////////////////////////////////////////////////////////////////
439 /// Return the current value:
440 /// the value of the PDF associated with the current index category state
441 
443 {
444  // Retrieve the proxy by index name
446 
447  //assert(proxy!=0) ;
448  if (proxy==0) return 0 ;
449 
450  // Calculate relative weighting factor for sim-pdfs of all extendable components
451  Double_t catFrac(1) ;
452  if (canBeExtended()) {
453  Double_t nEvtCat = ((RooAbsPdf*)(proxy->absArg()))->expectedEvents(_normSet) ;
454 
455  Double_t nEvtTot(0) ;
457  RooRealProxy* proxy2 ;
458  while((proxy2=(RooRealProxy*)iter->Next())) {
459  nEvtTot += ((RooAbsPdf*)(proxy2->absArg()))->expectedEvents(_normSet) ;
460  }
461  delete iter ;
462  catFrac=nEvtCat/nEvtTot ;
463  }
464 
465  // Return the selected PDF value, normalized by the number of index states
466  return ((RooAbsPdf*)(proxy->absArg()))->getVal(_normSet)*catFrac ;
467 }
468 
469 
470 
471 ////////////////////////////////////////////////////////////////////////////////
472 /// Return the number of expected events: If the index is in nset,
473 /// then return the sum of the expected events of all components,
474 /// otherwise return the number of expected events of the PDF
475 /// associated with the current index category state
476 
478 {
479  if (nset->contains(_indexCat.arg())) {
480 
481  Double_t sum(0) ;
482 
484  RooRealProxy* proxy ;
485  while((proxy=(RooRealProxy*)iter->Next())) {
486  sum += ((RooAbsPdf*)(proxy->absArg()))->expectedEvents(nset) ;
487  }
488  delete iter ;
489 
490  return sum ;
491 
492  } else {
493 
494  // Retrieve the proxy by index name
496 
497  //assert(proxy!=0) ;
498  if (proxy==0) return 0 ;
499 
500  // Return the selected PDF value, normalized by the number of index states
501  return ((RooAbsPdf*)(proxy->absArg()))->expectedEvents(nset);
502  }
503 }
504 
505 
506 
507 ////////////////////////////////////////////////////////////////////////////////
508 /// Forward determination of analytical integration capabilities to component p.d.f.s
509 /// A unique code is assigned to the combined integration capabilities of all associated
510 /// p.d.f.s
511 
513  const RooArgSet* normSet, const char* rangeName) const
514 {
515  // Declare that we can analytically integrate all requested observables
516  analVars.add(allVars) ;
517 
518  // Retrieve (or create) the required partial integral list
519  Int_t code ;
520 
521  // Check if this configuration was created before
522  CacheElem* cache = (CacheElem*) _partIntMgr.getObj(normSet,&analVars,0,RooNameReg::ptr(rangeName)) ;
523  if (cache) {
524  code = _partIntMgr.lastIndex() ;
525  return code+1 ;
526  }
527  cache = new CacheElem ;
528 
529  // Create the partial integral set for this request
531  RooRealProxy* proxy ;
532  while((proxy=(RooRealProxy*)iter->Next())) {
533  RooAbsReal* pdfInt = proxy->arg().createIntegral(analVars,normSet,0,rangeName) ;
534  cache->_partIntList.addOwned(*pdfInt) ;
535  }
536  delete iter ;
537 
538  // Store the partial integral list and return the assigned code ;
539  code = _partIntMgr.setObj(normSet,&analVars,cache,RooNameReg::ptr(rangeName)) ;
540 
541  return code+1 ;
542 }
543 
544 
545 
546 ////////////////////////////////////////////////////////////////////////////////
547 /// Return analytical integration defined by given code
548 
549 Double_t RooSimultaneous::analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* /*rangeName*/) const
550 {
551  // No integration scenario
552  if (code==0) {
553  return getVal(normSet) ;
554  }
555 
556  // Partial integration scenarios, rangeName already encoded in 'code'
557  CacheElem* cache = (CacheElem*) _partIntMgr.getObjByIndex(code-1) ;
558 
560  Int_t idx = _pdfProxyList.IndexOf(proxy) ;
561  return ((RooAbsReal*)cache->_partIntList.at(idx))->getVal(normSet) ;
562 }
563 
564 
565 
566 
567 
568 
569 ////////////////////////////////////////////////////////////////////////////////
570 /// Back-end for plotOn() implementation on RooSimultaneous which
571 /// needs special handling because a RooSimultaneous PDF cannot
572 /// project out its index category via integration, plotOn() will
573 /// abort if this is requested without providing a projection dataset
574 
576 {
577  // Sanity checks
578  if (plotSanityChecks(frame)) return frame ;
579 
580  // Extract projection configuration from command list
581  RooCmdConfig pc(Form("RooSimultaneous::plotOn(%s)",GetName())) ;
582  pc.defineString("sliceCatState","SliceCat",0,"",kTRUE) ;
583  pc.defineDouble("scaleFactor","Normalization",0,1.0) ;
584  pc.defineInt("scaleType","Normalization",0,RooAbsPdf::Relative) ;
585  pc.defineObject("sliceCatList","SliceCat",0,0,kTRUE) ;
586  pc.defineObject("projSet","Project",0) ;
587  pc.defineObject("sliceSet","SliceVars",0) ;
588  pc.defineObject("projDataSet","ProjData",0) ;
589  pc.defineObject("projData","ProjData",1) ;
590  pc.defineMutex("Project","SliceVars") ;
591  pc.allowUndefined() ; // there may be commands we don't handle here
592 
593  // Process and check varargs
594  pc.process(cmdList) ;
595  if (!pc.ok(kTRUE)) {
596  return frame ;
597  }
598 
599  const RooAbsData* projData = (const RooAbsData*) pc.getObject("projData") ;
600  const RooArgSet* projDataSet = (const RooArgSet*) pc.getObject("projDataSet") ;
601  const RooArgSet* sliceSetTmp = (const RooArgSet*) pc.getObject("sliceSet") ;
602  RooArgSet* sliceSet = sliceSetTmp ? ((RooArgSet*) sliceSetTmp->Clone()) : 0 ;
603  const RooArgSet* projSet = (const RooArgSet*) pc.getObject("projSet") ;
604  Double_t scaleFactor = pc.getDouble("scaleFactor") ;
605  ScaleType stype = (ScaleType) pc.getInt("scaleType") ;
606 
607 
608  // Look for category slice arguments and add them to the master slice list if found
609  const char* sliceCatState = pc.getString("sliceCatState",0,kTRUE) ;
610  const RooLinkedList& sliceCatList = pc.getObjectList("sliceCatList") ;
611  if (sliceCatState) {
612 
613  // Make the master slice set if it doesnt exist
614  if (!sliceSet) {
615  sliceSet = new RooArgSet ;
616  }
617 
618  // Prepare comma separated label list for parsing
619  char buf[1024] ;
620  strlcpy(buf,sliceCatState,1024) ;
621  const char* slabel = strtok(buf,",") ;
622 
623  // Loop over all categories provided by (multiple) Slice() arguments
624  TIterator* iter = sliceCatList.MakeIterator() ;
625  RooCategory* scat ;
626  while((scat=(RooCategory*)iter->Next())) {
627  if (slabel) {
628  // Set the slice position to the value indicate by slabel
629  scat->setLabel(slabel) ;
630  // Add the slice category to the master slice set
631  sliceSet->add(*scat,kFALSE) ;
632  }
633  slabel = strtok(0,",") ;
634  }
635  delete iter ;
636  }
637 
638  // Check if we have a projection dataset
639  if (!projData) {
640  coutE(InputArguments) << "RooSimultaneous::plotOn(" << GetName() << ") ERROR: must have a projection dataset for index category" << endl ;
641  return frame ;
642  }
643 
644  // Make list of variables to be projected
645  RooArgSet projectedVars ;
646  if (sliceSet) {
647  //cout << "frame->getNormVars() = " ; frame->getNormVars()->Print("1") ;
648 
649  makeProjectionSet(frame->getPlotVar(),frame->getNormVars(),projectedVars,kTRUE) ;
650 
651  // Take out the sliced variables
652  TIterator* iter = sliceSet->createIterator() ;
653  RooAbsArg* sliceArg ;
654  while((sliceArg=(RooAbsArg*)iter->Next())) {
655  RooAbsArg* arg = projectedVars.find(sliceArg->GetName()) ;
656  if (arg) {
657  projectedVars.remove(*arg) ;
658  } else {
659  coutI(Plotting) << "RooAbsReal::plotOn(" << GetName() << ") slice variable "
660  << sliceArg->GetName() << " was not projected anyway" << endl ;
661  }
662  }
663  delete iter ;
664  } else if (projSet) {
665  makeProjectionSet(frame->getPlotVar(),projSet,projectedVars,kFALSE) ;
666  } else {
667  makeProjectionSet(frame->getPlotVar(),frame->getNormVars(),projectedVars,kTRUE) ;
668  }
669 
670  Bool_t projIndex(kFALSE) ;
671 
672  if (!_indexCat.arg().isDerived()) {
673  // *** Error checking for a fundamental index category ***
674  //cout << "RooSim::plotOn: index is fundamental" << endl ;
675 
676  // Check that the provided projection dataset contains our index variable
677  if (!projData->get()->find(_indexCat.arg().GetName())) {
678  coutE(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") ERROR: Projection over index category "
679  << "requested, but projection data set doesn't contain index category" << endl ;
680  return frame ;
681  }
682 
683  if (projectedVars.find(_indexCat.arg().GetName())) {
684  projIndex=kTRUE ;
685  }
686 
687  } else {
688  // *** Error checking for a composite index category ***
689 
690  // Determine if any servers of the index category are in the projectedVars
691  TIterator* sIter = _indexCat.arg().serverIterator() ;
692  RooAbsArg* server ;
693  RooArgSet projIdxServers ;
694  Bool_t anyServers(kFALSE) ;
695  while((server=(RooAbsArg*)sIter->Next())) {
696  if (projectedVars.find(server->GetName())) {
697  anyServers=kTRUE ;
698  projIdxServers.add(*server) ;
699  }
700  }
701  delete sIter ;
702 
703  // Check that the projection dataset contains all the
704  // index category components we're projecting over
705 
706  // Determine if all projected servers of the index category are in the projection dataset
707  sIter = projIdxServers.createIterator() ;
708  Bool_t allServers(kTRUE) ;
709  while((server=(RooAbsArg*)sIter->Next())) {
710  if (!projData->get()->find(server->GetName())) {
711  allServers=kFALSE ;
712  }
713  }
714  delete sIter ;
715 
716  if (!allServers) {
717  coutE(Plotting) << "RooSimultaneous::plotOn(" << GetName()
718  << ") ERROR: Projection dataset doesn't contain complete set of index category dependents" << endl ;
719  return frame ;
720  }
721 
722  if (anyServers) {
723  projIndex = kTRUE ;
724  }
725  }
726 
727  // Calculate relative weight fractions of components
728  Roo1DTable* wTable = projData->table(_indexCat.arg()) ;
729 
730  // If we don't project over the index, just do the regular plotOn
731  if (!projIndex) {
732 
733  coutI(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") plot on " << frame->getPlotVar()->GetName()
734  << " represents a slice in the index category (" << _indexCat.arg().GetName() << ")" << endl ;
735 
736  // Reduce projData: take out fitCat (component) columns and entries that don't match selected slice
737  // Construct cut string to only select projection data event that match the current slice
738 
739  const RooAbsData* projDataTmp(projData) ;
740  if (projData) {
741  // Make list of categories columns to exclude from projection data
742  RooArgSet* indexCatComps = _indexCat.arg().getObservables(frame->getNormVars());
743 
744  // Make cut string to exclude rows from projection data
745  TString cutString ;
746  TIterator* compIter = indexCatComps->createIterator() ;
747  RooAbsCategory* idxComp ;
748  Bool_t first(kTRUE) ;
749  while((idxComp=(RooAbsCategory*)compIter->Next())) {
750  if (!first) {
751  cutString.Append("&&") ;
752  } else {
753  first=kFALSE ;
754  }
755  cutString.Append(Form("%s==%d",idxComp->GetName(),idxComp->getIndex())) ;
756  }
757  delete compIter ;
758 
759  // Make temporary projData without RooSim index category components
760  RooArgSet projDataVars(*projData->get()) ;
761  projDataVars.remove(*indexCatComps,kTRUE,kTRUE) ;
762 
763  projDataTmp = ((RooAbsData*)projData)->reduce(projDataVars,cutString) ;
764  delete indexCatComps ;
765  }
766 
767  // Multiply scale factor with fraction of events in current state of index
768 
769 // RooPlot* retFrame = getPdf(_indexCat.arg().getLabel())->plotOn(frame,drawOptions,
770 // scaleFactor*wTable->getFrac(_indexCat.arg().getLabel()),
771 // stype,projDataTmp,projSet) ;
772 
773  // Override normalization and projection dataset
774  RooCmdArg tmp1 = RooFit::Normalization(scaleFactor*wTable->getFrac(_indexCat.arg().getLabel()),stype) ;
775  RooCmdArg tmp2 = RooFit::ProjWData(*projDataSet,*projDataTmp) ;
776 
777  // WVE -- do not adjust normalization for asymmetry plots
778  RooLinkedList cmdList2(cmdList) ;
779  if (!cmdList.find("Asymmetry")) {
780  cmdList2.Add(&tmp1) ;
781  }
782  cmdList2.Add(&tmp2) ;
783 
784  // Plot single component
785  RooPlot* retFrame = getPdf(_indexCat.arg().getLabel())->plotOn(frame,cmdList2) ;
786 
787  // Delete temporary dataset
788  if (projDataTmp) {
789  delete projDataTmp ;
790  }
791 
792  delete wTable ;
793  delete sliceSet ;
794  return retFrame ;
795  }
796 
797  // If we project over the index, plot using a temporary RooAddPdf
798  // using the weights from the data as coefficients
799 
800  // Make a deep clone of our index category
801  RooArgSet* idxCloneSet = (RooArgSet*) RooArgSet(_indexCat.arg()).snapshot(kTRUE) ;
802  RooAbsCategoryLValue* idxCatClone = (RooAbsCategoryLValue*) idxCloneSet->find(_indexCat.arg().GetName()) ;
803 
804  // Build the list of indexCat components that are sliced
805  RooArgSet* idxCompSliceSet = _indexCat.arg().getObservables(frame->getNormVars()) ;
806  idxCompSliceSet->remove(projectedVars,kTRUE,kTRUE) ;
807  TIterator* idxCompSliceIter = idxCompSliceSet->createIterator() ;
808 
809  // Make a new expression that is the weighted sum of requested components
810  RooArgList pdfCompList ;
811  RooArgList wgtCompList ;
812 //RooAbsPdf* pdf ;
813  RooRealProxy* proxy ;
815  Double_t sumWeight(0) ;
816  while((proxy=(RooRealProxy*)pIter->Next())) {
817 
818  idxCatClone->setLabel(proxy->name()) ;
819 
820  // Determine if this component is the current slice (if we slice)
821  Bool_t skip(kFALSE) ;
822  idxCompSliceIter->Reset() ;
823  RooAbsCategory* idxSliceComp ;
824  while((idxSliceComp=(RooAbsCategory*)idxCompSliceIter->Next())) {
825  RooAbsCategory* idxComp = (RooAbsCategory*) idxCloneSet->find(idxSliceComp->GetName()) ;
826  if (idxComp->getIndex()!=idxSliceComp->getIndex()) {
827  skip=kTRUE ;
828  break ;
829  }
830  }
831  if (skip) continue ;
832 
833  // Instantiate a RRV holding this pdfs weight fraction
834  RooRealVar *wgtVar = new RooRealVar(proxy->name(),"coef",wTable->getFrac(proxy->name())) ;
835  wgtCompList.addOwned(*wgtVar) ;
836  sumWeight += wTable->getFrac(proxy->name()) ;
837 
838  // Add the PDF to list list
839  pdfCompList.add(proxy->arg()) ;
840  }
841 
842  TString plotVarName(GetName()) ;
843  RooAddPdf *plotVar = new RooAddPdf(plotVarName,"weighted sum of RS components",pdfCompList,wgtCompList) ;
844 
845  // Fix appropriate coefficient normalization in plot function
846  if (_plotCoefNormSet.getSize()>0) {
848  }
849 
850  RooAbsData* projDataTmp(0) ;
851  RooArgSet projSetTmp ;
852  if (projData) {
853 
854  // Construct cut string to only select projection data event that match the current slice
855  TString cutString ;
856  if (idxCompSliceSet->getSize()>0) {
857  idxCompSliceIter->Reset() ;
858  RooAbsCategory* idxSliceComp ;
859  Bool_t first(kTRUE) ;
860  while((idxSliceComp=(RooAbsCategory*)idxCompSliceIter->Next())) {
861  if (!first) {
862  cutString.Append("&&") ;
863  } else {
864  first=kFALSE ;
865  }
866  cutString.Append(Form("%s==%d",idxSliceComp->GetName(),idxSliceComp->getIndex())) ;
867  }
868  }
869 
870  // Make temporary projData without RooSim index category components
871  RooArgSet projDataVars(*projData->get()) ;
872  RooArgSet* idxCatServers = _indexCat.arg().getObservables(frame->getNormVars()) ;
873 
874  projDataVars.remove(*idxCatServers,kTRUE,kTRUE) ;
875 
876  if (idxCompSliceSet->getSize()>0) {
877  projDataTmp = ((RooAbsData*)projData)->reduce(projDataVars,cutString) ;
878  } else {
879  projDataTmp = ((RooAbsData*)projData)->reduce(projDataVars) ;
880  }
881 
882 
883 
884  if (projSet) {
885  projSetTmp.add(*projSet) ;
886  projSetTmp.remove(*idxCatServers,kTRUE,kTRUE);
887  }
888 
889 
890  delete idxCatServers ;
891  }
892 
893 
894  if (_indexCat.arg().isDerived() && idxCompSliceSet->getSize()>0) {
895  coutI(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") plot on " << frame->getPlotVar()->GetName()
896  << " represents a slice in index category components " << *idxCompSliceSet << endl ;
897 
898  RooArgSet* idxCompProjSet = _indexCat.arg().getObservables(frame->getNormVars()) ;
899  idxCompProjSet->remove(*idxCompSliceSet,kTRUE,kTRUE) ;
900  if (idxCompProjSet->getSize()>0) {
901  coutI(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") plot on " << frame->getPlotVar()->GetName()
902  << " averages with data index category components " << *idxCompProjSet << endl ;
903  }
904  delete idxCompProjSet ;
905  } else {
906  coutI(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") plot on " << frame->getPlotVar()->GetName()
907  << " averages with data index category (" << _indexCat.arg().GetName() << ")" << endl ;
908  }
909 
910 
911  // Override normalization and projection dataset
912  RooLinkedList cmdList2(cmdList) ;
913 
914  RooCmdArg tmp1 = RooFit::Normalization(scaleFactor*sumWeight,stype) ;
915  RooCmdArg tmp2 = RooFit::ProjWData(*projDataSet,*projDataTmp) ;
916  // WVE -- do not adjust normalization for asymmetry plots
917  if (!cmdList.find("Asymmetry")) {
918  cmdList2.Add(&tmp1) ;
919  }
920  cmdList2.Add(&tmp2) ;
921 
922  RooPlot* frame2 ;
923  if (projSetTmp.getSize()>0) {
924  // Plot temporary function
925  RooCmdArg tmp3 = RooFit::Project(projSetTmp) ;
926  cmdList2.Add(&tmp3) ;
927  frame2 = plotVar->plotOn(frame,cmdList2) ;
928  } else {
929  // Plot temporary function
930  frame2 = plotVar->plotOn(frame,cmdList2) ;
931  }
932 
933  // Cleanup
934  delete sliceSet ;
935  delete pIter ;
936  delete wTable ;
937  delete idxCloneSet ;
938  delete idxCompSliceIter ;
939  delete idxCompSliceSet ;
940  delete plotVar ;
941 
942  if (projDataTmp) delete projDataTmp ;
943 
944  return frame2 ;
945 }
946 
947 
948 
949 ////////////////////////////////////////////////////////////////////////////////
950 /// OBSOLETE -- Retained for backward compatibility
951 
952 RooPlot* RooSimultaneous::plotOn(RooPlot *frame, Option_t* drawOptions, Double_t scaleFactor,
953  ScaleType stype, const RooAbsData* projData, const RooArgSet* projSet,
954  Double_t /*precision*/, Bool_t /*shiftToZero*/, const RooArgSet* /*projDataSet*/,
955  Double_t /*rangeLo*/, Double_t /*rangeHi*/, RooCurve::WingMode /*wmode*/) const
956 {
957  // Make command list
958  RooLinkedList cmdList ;
959  cmdList.Add(new RooCmdArg(RooFit::DrawOption(drawOptions))) ;
960  cmdList.Add(new RooCmdArg(RooFit::Normalization(scaleFactor,stype))) ;
961  if (projData) cmdList.Add(new RooCmdArg(RooFit::ProjWData(*projData))) ;
962  if (projSet) cmdList.Add(new RooCmdArg(RooFit::Project(*projSet))) ;
963 
964  // Call new method
965  RooPlot* ret = plotOn(frame,cmdList) ;
966 
967  // Cleanup
968  cmdList.Delete() ;
969  return ret ;
970 }
971 
972 
973 
974 ////////////////////////////////////////////////////////////////////////////////
975 /// Interface function used by test statistics to freeze choice of observables
976 /// for interpretation of fraction coefficients. Needed here because a RooSimultaneous
977 /// works like a RooAddPdf when plotted
978 
979 void RooSimultaneous::selectNormalization(const RooArgSet* normSet, Bool_t /*force*/)
980 {
982  if (normSet) _plotCoefNormSet.add(*normSet) ;
983 }
984 
985 
986 ////////////////////////////////////////////////////////////////////////////////
987 /// Interface function used by test statistics to freeze choice of range
988 /// for interpretation of fraction coefficients. Needed here because a RooSimultaneous
989 /// works like a RooAddPdf when plotted
990 
991 void RooSimultaneous::selectNormalizationRange(const char* normRange2, Bool_t /*force*/)
992 {
993  _plotCoefNormRange = RooNameReg::ptr(normRange2) ;
994 }
995 
996 
997 
998 
999 ////////////////////////////////////////////////////////////////////////////////
1000 
1002  const RooArgSet* auxProto, Bool_t verbose, Bool_t autoBinned, const char* binnedTag) const
1003 {
1004  const char* idxCatName = _indexCat.arg().GetName() ;
1005 
1006  if (vars.find(idxCatName) && prototype==0 && (auxProto==0 || auxProto->getSize()==0) && (autoBinned || (binnedTag && strlen(binnedTag)))) {
1007 
1008  // Return special generator config that can also do binned generation for selected states
1009  return new RooSimSplitGenContext(*this,vars,verbose,autoBinned,binnedTag) ;
1010 
1011  } else {
1012 
1013  // Return regular generator config ;
1014  return genContext(vars,prototype,auxProto,verbose) ;
1015  }
1016 }
1017 
1018 
1019 
1020 ////////////////////////////////////////////////////////////////////////////////
1021 /// Return specialized generator contenxt for simultaneous p.d.f.s
1022 
1024  const RooArgSet* auxProto, Bool_t verbose) const
1025 {
1026  const char* idxCatName = _indexCat.arg().GetName() ;
1027  const RooArgSet* protoVars = prototype ? prototype->get() : 0 ;
1028 
1029  if (vars.find(idxCatName) || (protoVars && protoVars->find(idxCatName))) {
1030 
1031  // Generating index category: return special sim-context
1032  return new RooSimGenContext(*this,vars,prototype,auxProto,verbose) ;
1033 
1034  } else if (_indexCat.arg().isDerived()) {
1035  // Generating dependents of a derived index category
1036 
1037  // Determine if we none,any or all servers
1038  Bool_t anyServer(kFALSE), allServers(kTRUE) ;
1039  if (prototype) {
1040  TIterator* sIter = _indexCat.arg().serverIterator() ;
1041  RooAbsArg* server ;
1042  while((server=(RooAbsArg*)sIter->Next())) {
1043  if (prototype->get()->find(server->GetName())) {
1044  anyServer=kTRUE ;
1045  } else {
1046  allServers=kFALSE ;
1047  }
1048  }
1049  delete sIter ;
1050  } else {
1051  allServers=kTRUE ;
1052  }
1053 
1054  if (allServers) {
1055  // Use simcontext if we have all servers
1056 
1057  return new RooSimGenContext(*this,vars,prototype,auxProto,verbose) ;
1058  } else if (!allServers && anyServer) {
1059  // Abort if we have only part of the servers
1060  coutE(Plotting) << "RooSimultaneous::genContext: ERROR: prototype must include either all "
1061  << " components of the RooSimultaneous index category or none " << endl ;
1062  return 0 ;
1063  }
1064  // Otherwise make single gencontext for current state
1065  }
1066 
1067  // Not generating index cat: return context for pdf associated with present index state
1069  if (!proxy) {
1070  coutE(InputArguments) << "RooSimultaneous::genContext(" << GetName()
1071  << ") ERROR: no PDF associated with current state ("
1072  << _indexCat.arg().GetName() << "=" << _indexCat.arg().getLabel() << ")" << endl ;
1073  return 0 ;
1074  }
1075  return ((RooAbsPdf*)proxy->absArg())->genContext(vars,prototype,auxProto,verbose) ;
1076 }
1077 
1078 
1079 
1080 
1081 ////////////////////////////////////////////////////////////////////////////////
1082 
1084  const RooArgSet* nset,
1085  Double_t scaleFactor,
1086  Bool_t correctForBinVolume,
1087  Bool_t showProgress) const
1088 {
1089  if (RooAbsReal::fillDataHist (hist, nset, scaleFactor,
1090  correctForBinVolume, showProgress) == 0)
1091  return 0;
1092 
1093  Double_t sum = 0;
1094  for (int i=0 ; i<hist->numEntries() ; i++) {
1095  hist->get(i) ;
1096  sum += hist->weight();
1097  }
1098  if (sum != 0) {
1099  for (int i=0 ; i<hist->numEntries() ; i++) {
1100  hist->get(i) ;
1101  hist->set (hist->weight() / sum);
1102  }
1103  }
1104 
1105  return hist;
1106 }
1107 
1108 
1109 
1110 
1111 ////////////////////////////////////////////////////////////////////////////////
1112 /// Special generator interface for generation of 'global observables' -- for RooStats tools
1113 
1115 {
1116  // Make set with clone of variables (placeholder for output)
1117  RooArgSet* globClone = (RooArgSet*) whatVars.snapshot() ;
1118 
1119  RooDataSet* data = new RooDataSet("gensimglobal","gensimglobal",whatVars) ;
1120 
1121  // Construct iterator over index types
1122  TIterator* iter = indexCat().typeIterator() ;
1123 
1124  for (Int_t i=0 ; i<nEvents ; i++) {
1125  iter->Reset() ;
1126  RooCatType* tt ;
1127  while((tt=(RooCatType*) iter->Next())) {
1128 
1129  // Get pdf associated with state from simpdf
1130  RooAbsPdf* pdftmp = getPdf(tt->GetName()) ;
1131 
1132  // Generate only global variables defined by the pdf associated with this state
1133  RooArgSet* globtmp = pdftmp->getObservables(whatVars) ;
1134  RooDataSet* tmp = pdftmp->generate(*globtmp,1) ;
1135 
1136  // Transfer values to output placeholder
1137  *globClone = *tmp->get(0) ;
1138 
1139  // Cleanup
1140  delete globtmp ;
1141  delete tmp ;
1142  }
1143  data->add(*globClone) ;
1144  }
1145 
1146 
1147  delete iter ;
1148  delete globClone ;
1149  return data ;
1150 }
1151 
1152 
1153 
1154 
1155 
1156 
1157 
1158 
1159 
RooArgSet * getVariables(Bool_t stripDisconnected=kTRUE) const
Return RooArgSet with all variables (tree leaf nodes of expresssion tree)
Definition: RooAbsArg.cxx:2082
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
const char * label() const
RooAddPdf is an efficient implementation of a sum of PDFs of the form.
Definition: RooAddPdf.h:29
TIterator * createIterator(Bool_t dir=kIterForward) const
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 fraction c...
RooSetProxy _plotCoefNormSet
static long int sum(long int i)
Definition: Factory.cxx:2162
#define coutE(a)
Definition: RooMsgService.h:34
virtual void selectNormalization(const RooArgSet *depSet=0, Bool_t force=kFALSE)
Interface function used by test statistics to freeze choice of observables for interpretation of frac...
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
virtual void Delete(Option_t *option="")
Remove all objects from the list AND delete all heap based objects.
Definition: TList.cxx:409
RooCmdArg DrawOption(const char *opt)
RooCategoryProxy _indexCat
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
virtual const RooArgSet * get() const
Definition: RooDataHist.h:77
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
virtual void Reset()=0
Bool_t defineDouble(const char *name, const char *argName, Int_t doubleNum, Double_t defValue=0.)
Define Double_t property name &#39;name&#39; mapped to Double_t in slot &#39;doubleNum&#39; in RooCmdArg with name ar...
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Definition: RooAbsArg.h:194
virtual const RooArgSet * get() const
Definition: RooAbsData.h:77
const char Option_t
Definition: RtypesCore.h:62
virtual Bool_t addOwned(RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
#define coutI(a)
Definition: RooMsgService.h:31
const char * getString(const char *name, const char *defaultValue="", Bool_t convEmptyToNull=kFALSE)
Return string property registered with name &#39;name&#39;.
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:79
virtual RooAbsGenContext * genContext(const RooArgSet &vars, const RooDataSet *prototype=0, const RooArgSet *auxProto=0, Bool_t verbose=kFALSE) const
Return specialized generator contenxt for simultaneous p.d.f.s.
TObject * find(const char *name) const
Return pointer to object with given name in collection.
void set(Double_t weight, Double_t wgtErr=-1)
Increment the weight of the bin enclosing the coordinates given by &#39;row&#39; by the specified amount...
virtual Bool_t setLabel(const char *label, Bool_t printError=kTRUE)=0
Double_t analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=0) const
Return analytical integration defined by given code.
T * getObjByIndex(Int_t index) const
const TNamed * _plotCoefNormRange
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
RooCmdArg Project(const RooArgSet &projSet)
RooAbsData * reduce(const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg(), const RooCmdArg &arg3=RooCmdArg(), const RooCmdArg &arg4=RooCmdArg(), const RooCmdArg &arg5=RooCmdArg(), const RooCmdArg &arg6=RooCmdArg(), const RooCmdArg &arg7=RooCmdArg(), const RooCmdArg &arg8=RooCmdArg())
Create a reduced copy of this dataset.
Definition: RooAbsData.cxx:343
virtual ExtendMode extendMode() const
WVE NEEDS FIX.
RooCmdArg ProjWData(const RooAbsData &projData, Bool_t binData=kFALSE)
virtual Double_t evaluate() const
Return the current value: the value of the PDF associated with the current index category state...
void allowUndefined(Bool_t flag=kTRUE)
Definition: RooCmdConfig.h:39
virtual Int_t getIndex() const
Return index number of current state.
virtual RooDataHist * fillDataHist(RooDataHist *hist, const RooArgSet *nset, Double_t scaleFactor, Bool_t correctForBinVolume=kFALSE, Bool_t showProgress=kFALSE) const
Bool_t addOwnedComponents(const RooArgSet &comps)
Take ownership of the contents of &#39;comps&#39;.
Definition: RooAbsArg.cxx:2282
STL namespace.
virtual RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1=RooCmdArg::none(), 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 RooCmdArg &arg9=RooCmdArg::none(), const RooCmdArg &arg10=RooCmdArg::none()) const
Plot (project) PDF on specified frame.
Definition: RooAbsPdf.h:105
#define coutW(a)
Definition: RooMsgService.h:33
TObject * getObject(const char *name, TObject *obj=0)
Return TObject property registered with name &#39;name&#39;.
Iterator abstract base class.
Definition: TIterator.h:30
Bool_t mustBeExtended() const
Definition: RooAbsPdf.h:220
virtual TObject * FindObject(const char *name) const
Find an object in this list using its name.
Definition: TList.cxx:501
virtual Bool_t setLabel(const char *label, Bool_t printError=kTRUE)
Set value by specifying the name of the desired state If printError is set, a message will be printed...
RooObjCacheManager _partIntMgr
RooDataHist * fillDataHist(RooDataHist *hist, const RooArgSet *nset, Double_t scaleFactor, Bool_t correctForBinVolume=kFALSE, Bool_t showProgress=kFALSE) const
Fill a RooDataHist with values sampled from this function at the bin centers.
Bool_t process(const RooCmdArg &arg)
Process given RooCmdArg.
RooDataSet is a container class to hold N-dimensional binned data.
Definition: RooDataHist.h:40
friend class CacheElem
Definition: RooAbsPdf.h:314
Int_t numTypes(const char *=0) const
virtual TObject * Clone(const char *newname=0) const
Make a clone of an object using the Streamer facility.
TText * tt
Definition: textangle.C:16
Bool_t addPdf(const RooAbsPdf &pdf, const char *catLabel)
Associate given PDF with index category state label &#39;catLabel&#39;.
void initialize(RooAbsCategoryLValue &inIndexCat, std::map< std::string, RooAbsPdf *> pdfMap)
T * getObj(const RooArgSet *nset, Int_t *sterileIndex=0, const TNamed *isetRangeName=0)
Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &numVars, const RooArgSet *normSet, const char *rangeName=0) const
Forward determination of analytical integration capabilities to component p.d.f.s A unique code is as...
Bool_t defineString(const char *name, const char *argName, Int_t stringNum, const char *defValue="", Bool_t appendMode=kFALSE)
Define Double_t property name &#39;name&#39; mapped to Double_t in slot &#39;stringNum&#39; in RooCmdArg with name ar...
const RooAbsCategory & arg() const
virtual Double_t weight() const
Definition: RooDataHist.h:96
RooCatType is an auxilary class for RooAbsCategory and defines a a single category state...
Definition: RooCatType.h:22
friend class RooArgSet
Definition: RooAbsArg.h:469
virtual RooAbsGenContext * autoGenContext(const RooArgSet &vars, const RooDataSet *prototype=0, const RooArgSet *auxProto=0, Bool_t verbose=kFALSE, Bool_t autoBinned=kTRUE, const char *binnedTag="") const
const RooArgSet * getNormVars() const
Definition: RooPlot.h:139
RooAbsGenContext is the abstract base class for generator contexts of RooAbsPdf objects.
RooAbsCategoryLValue is the common abstract base class for objects that represent a discrete value th...
Bool_t defineInt(const char *name, const char *argName, Int_t intNum, Int_t defValue=0)
Define integer property name &#39;name&#39; mapped to integer in slot &#39;intNum&#39; in RooCmdArg with name argName...
virtual TIterator * MakeIterator(Bool_t dir=kIterForward) const
Return a list iterator.
Definition: TList.cxx:608
virtual Double_t expectedEvents(const RooArgSet *nset) const
Return the number of expected events: If the index is in nset, then return the sum of the expected ev...
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:36
virtual Bool_t isFundamental() const
Definition: RooAbsArg.h:157
virtual Bool_t setArg(RooAbsCategory &newRef)
Change object held in proxy into newRef.
virtual const char * name() const
Definition: RooArgProxy.h:42
virtual void Add(TObject *arg)
Definition: RooLinkedList.h:62
const int nEvents
Definition: testRooFit.cxx:42
void makeProjectionSet(const RooAbsArg *plotVar, const RooArgSet *allVars, RooArgSet &projectedVars, Bool_t silent) const
Utility function for plotOn() that constructs the set of observables to project when plotting ourselv...
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:501
void defineMutex(const char *argName1, const char *argName2)
Define arguments named argName1 and argName2 mutually exclusive.
Int_t getInt(const char *name, Int_t defaultValue=0)
Return integer property registered with name &#39;name&#39;.
RooAbsCollection * snapshot(Bool_t deepCopy=kTRUE) const
Take a snap shot of current collection contents: An owning collection is returned containing clones o...
friend class RooSimSplitGenContext
TIterator * typeIterator() const
Return iterator over all defined states.
RooAbsArg * at(Int_t idx) const
Definition: RooArgList.h:84
virtual const Text_t * GetName() const
Returns name of object.
Definition: RooCatType.h:44
const RooAbsCategoryLValue & indexCat() const
static const TNamed * ptr(const char *stringPtr)
Return a unique TNamed pointer for given C++ string.
Definition: RooNameReg.cxx:125
Bool_t ok(Bool_t verbose) const
Return true of parsing was successful.
bool verbose
char * Form(const char *fmt,...)
Int_t lastIndex() const
RooAbsArg * absArg() const
Definition: RooArgProxy.h:37
Class RooCmdConfig is a configurable parser for RooCmdArg named arguments.
Definition: RooCmdConfig.h:27
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition: RooAbsData.h:37
virtual const char * getLabel() const
Return label string of current state.
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:29
RooCategory represents a fundamental (non-derived) discrete value object.
Definition: RooCategory.h:24
Bool_t canBeExtended() const
Definition: RooAbsPdf.h:216
TIterator * serverIterator() const
Definition: RooAbsArg.h:112
virtual const RooArgSet * get(Int_t index) const
Return RooArgSet with coordinates of event &#39;index&#39;.
virtual void add(const RooArgSet &row, Double_t weight=1.0, Double_t weightError=0)
Add a data point, with its coordinates specified in the &#39;data&#39; argset, to the data set...
A RooPlot is a plot frame and a container for graphics objects within that frame. ...
Definition: RooPlot.h:41
const Bool_t kFALSE
Definition: RtypesCore.h:92
RooLinkedList is an collection class for internal use, storing a collection of RooAbsArg pointers in ...
Definition: RooLinkedList.h:35
void Delete(Option_t *o=0)
Remove all elements in collection and delete all elements NB: Collection does not own elements...
friend class RooAddPdf
Definition: RooAbsReal.h:471
#define ClassImp(name)
Definition: Rtypes.h:336
RooAbsArg * find(const char *name) const
Find object with given name in list.
virtual Int_t numEntries() const
Return the number of bins.
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
friend class RooDataSet
Definition: RooAbsArg.h:534
Double_t getDouble(const char *name, Double_t defaultValue=0)
Return Double_t property registered with name &#39;name&#39;.
int type
Definition: TGX11.cxx:120
RooCmdArg Normalization(Double_t scaleFactor)
const RooLinkedList & getObjectList(const char *name)
Return list of objects registered with name &#39;name&#39;.
Mother of all ROOT objects.
Definition: TObject.h:37
virtual RooDataSet * generateSimGlobal(const RooArgSet &whatVars, Int_t nEvents)
Special generator interface for generation of &#39;global observables&#39; – for RooStats tools...
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE)
Remove the specified argument from our list.
RooAbsPdf * getPdf(const char *catName) const
Return the p.d.f associated with the given index category name.
virtual ~RooSimultaneous()
Destructor.
friend class RooRealProxy
Definition: RooAbsReal.h:403
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Definition: RooAbsPdf.h:41
virtual Bool_t isDerived() const
Definition: RooAbsArg.h:81
virtual void Add(TObject *obj)
Definition: TList.h:77
RooAbsRealLValue * getPlotVar() const
Definition: RooPlot.h:132
TIterator * MakeIterator(Bool_t dir=kTRUE) const
Return an iterator over this list.
RooDataSet * generate(const RooArgSet &whatVars, Int_t nEvents, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none())
Generate a new dataset containing the specified variables with events sampled from our distribution...
Definition: RooAbsPdf.cxx:1702
virtual TObject * Next()=0
RooAbsCategory is the common abstract base class for objects that represent a discrete value with a f...
virtual RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg(), const RooCmdArg &arg3=RooCmdArg(), const RooCmdArg &arg4=RooCmdArg(), const RooCmdArg &arg5=RooCmdArg(), const RooCmdArg &arg6=RooCmdArg(), const RooCmdArg &arg7=RooCmdArg(), const RooCmdArg &arg8=RooCmdArg(), const RooCmdArg &arg9=RooCmdArg(), const RooCmdArg &arg10=RooCmdArg()) const
Plot (project) PDF on specified frame.
RooRealProxy is the concrete proxy for RooAbsReal objects A RooRealProxy is the general mechanism to ...
Definition: RooRealProxy.h:23
RooSuperCategory consolidates several RooAbsCategoryLValue objects into a single category.
Bool_t defineObject(const char *name, const char *argName, Int_t setNum, const TObject *obj=0, Bool_t isArray=kFALSE)
Define TObject property name &#39;name&#39; mapped to object in slot &#39;setNum&#39; in RooCmdArg with name argName ...
RooArgSet * _normSet
Normalization integral (owned by _normMgr)
Definition: RooAbsPdf.h:302
Int_t setObj(const RooArgSet *nset, T *obj, const TNamed *isetRangeName=0)
Bool_t contains(const RooAbsArg &var) const
Definition: first.py:1
virtual void fixAddCoefNormalization(const RooArgSet &addNormSet=RooArgSet(), Bool_t force=kTRUE)
Fix the interpretation of the coefficient of any RooAddPdf component in the expression tree headed by...
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
virtual Bool_t setLabel(const char *label, Bool_t printError=kTRUE)
Set the value of the super category by specifying the state name by setting the state names of the co...
RooSimultaneous facilitates simultaneous fitting of multiple PDFs to subsets of a given dataset...
const Bool_t kTRUE
Definition: RtypesCore.h:91
Double_t getFrac(const char *label, Bool_t silent=kFALSE) const
Return the fraction of entries in the table contained in the slot named &#39;label&#39;.
Definition: Roo1DTable.cxx:301
friend class RooSimGenContext
static char * skip(char **buf, const char *delimiters)
Definition: civetweb.c:2039
Bool_t plotSanityChecks(RooPlot *frame) const
Utility function for plotOn(), perform general sanity check on frame to ensure safe plotting operatio...
Roo1DTable implements a one-dimensional table.
Definition: Roo1DTable.h:24
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...
RooCmdArg is a named container for two doubles, two integers two object points and three string point...
Definition: RooCmdArg.h:27
virtual Int_t IndexOf(const TObject *obj) const
Return index of object in collection.
virtual Roo1DTable * table(const RooArgSet &catSet, const char *cuts="", const char *opts="") const
Construct table for product of categories in catSet.
Definition: RooAbsData.cxx:752