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