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