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,bool,bool) 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 "RooSimultaneous.h"
51#include "RooPlot.h"
52#include "RooCurve.h"
53#include "RooRealVar.h"
54#include "RooAddPdf.h"
55#include "RooAbsData.h"
56#include "Roo1DTable.h"
57#include "RooSimGenContext.h"
59#include "RooDataSet.h"
60#include "RooCmdConfig.h"
61#include "RooNameReg.h"
62#include "RooGlobalFunc.h"
63#include "RooMsgService.h"
64#include "RooCategory.h"
65#include "RooSuperCategory.h"
66#include "RooDataHist.h"
67#include "RooRandom.h"
68#include "RooArgSet.h"
69#include "RooBinSamplingPdf.h"
70
71#include "ROOT/StringUtils.hxx"
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,false,false),
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,false,false),
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,false,false),
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 {
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 simComps(false) ;
172 for (map<string,RooAbsPdf*>::iterator iter=pdfMap.begin() ; iter!=pdfMap.end() ; ++iter) {
173 if (dynamic_cast<RooSimultaneous*>(iter->second)) {
174 simComps = true ;
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) {
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),true) ;
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.empty()) {
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 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 true ;
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 true ;
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 true ;
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 false ;
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 anyCanExtend(false) ;
420 bool anyMustExtend(false) ;
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 = true;
431 }
432 if (pdf->mustBeExtended()) {
433 //cout << "RooSim::extendedMode(" << GetName() << ") component " << pdf->GetName() << " MUST be extended" << endl;
434 anyMustExtend = true;
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 catFrac(1) ;
467 if (canBeExtended()) {
468 double nEvtCat = ((RooAbsPdf*)(proxy->absArg()))->expectedEvents(_normSet) ;
469
470 double 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 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 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,"",true) ;
598 pc.defineDouble("scaleFactor","Normalization",0,1.0) ;
599 pc.defineInt("scaleType","Normalization",0,RooAbsPdf::Relative) ;
600 pc.defineObject("sliceCatList","SliceCat",0,0,true) ;
601 // This dummy is needed for plotOn to recognize the "SliceCatMany" command.
602 // It is not used directly, but the "SliceCat" commands are nested in it.
603 // Removing this dummy definition results in "ERROR: unrecognized command: SliceCatMany".
604 pc.defineObject("dummy1","SliceCatMany",0) ;
605 pc.defineObject("projSet","Project",0) ;
606 pc.defineObject("sliceSet","SliceVars",0) ;
607 pc.defineObject("projDataSet","ProjData",0) ;
608 pc.defineObject("projData","ProjData",1) ;
609 pc.defineMutex("Project","SliceVars") ;
610 pc.allowUndefined() ; // there may be commands we don't handle here
611
612 // Process and check varargs
613 pc.process(cmdList) ;
614 if (!pc.ok(true)) {
615 return frame ;
616 }
617
618 const RooAbsData* projData = (const RooAbsData*) pc.getObject("projData") ;
619 const RooArgSet* projDataSet = (const RooArgSet*) pc.getObject("projDataSet") ;
620 const RooArgSet* sliceSetTmp = (const RooArgSet*) pc.getObject("sliceSet") ;
621 std::unique_ptr<RooArgSet> sliceSet( sliceSetTmp ? ((RooArgSet*) sliceSetTmp->Clone()) : nullptr );
622 const RooArgSet* projSet = (const RooArgSet*) pc.getObject("projSet") ;
623 double scaleFactor = pc.getDouble("scaleFactor") ;
624 ScaleType stype = (ScaleType) pc.getInt("scaleType") ;
625
626
627 // Look for category slice arguments and add them to the master slice list if found
628 const char* sliceCatState = pc.getString("sliceCatState",0,true) ;
629 const RooLinkedList& sliceCatList = pc.getObjectList("sliceCatList") ;
630 if (sliceCatState) {
631
632 // Make the master slice set if it doesnt exist
633 if (!sliceSet) {
634 sliceSet.reset(new RooArgSet);
635 }
636
637 // Prepare comma separated label list for parsing
638 auto catTokens = ROOT::Split(sliceCatState, ",");
639
640 // Loop over all categories provided by (multiple) Slice() arguments
641 unsigned int tokenIndex = 0;
642 for(auto * scat : static_range_cast<RooCategory*>(sliceCatList)) {
643 const char* slabel = tokenIndex >= catTokens.size() ? nullptr : catTokens[tokenIndex++].c_str();
644
645 if (slabel) {
646 // Set the slice position to the value indicated by slabel
647 scat->setLabel(slabel) ;
648 // Add the slice category to the master slice set
649 sliceSet->add(*scat,false) ;
650 }
651 }
652 }
653
654 // Check if we have a projection dataset
655 if (!projData) {
656 coutE(InputArguments) << "RooSimultaneous::plotOn(" << GetName() << ") ERROR: must have a projection dataset for index category" << endl ;
657 return frame ;
658 }
659
660 // Make list of variables to be projected
661 RooArgSet projectedVars ;
662 if (sliceSet) {
663 makeProjectionSet(frame->getPlotVar(),frame->getNormVars(),projectedVars,true) ;
664
665 // Take out the sliced variables
666 for (const auto sliceArg : *sliceSet) {
667 RooAbsArg* arg = projectedVars.find(sliceArg->GetName()) ;
668 if (arg) {
669 projectedVars.remove(*arg) ;
670 } else {
671 coutI(Plotting) << "RooAbsReal::plotOn(" << GetName() << ") slice variable "
672 << sliceArg->GetName() << " was not projected anyway" << endl ;
673 }
674 }
675 } else if (projSet) {
676 makeProjectionSet(frame->getPlotVar(),projSet,projectedVars,false) ;
677 } else {
678 makeProjectionSet(frame->getPlotVar(),frame->getNormVars(),projectedVars,true) ;
679 }
680
681 bool projIndex(false) ;
682
683 if (!_indexCat.arg().isDerived()) {
684 // *** Error checking for a fundamental index category ***
685 //cout << "RooSim::plotOn: index is fundamental" << endl ;
686
687 // Check that the provided projection dataset contains our index variable
688 if (!projData->get()->find(_indexCat.arg().GetName())) {
689 coutE(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") ERROR: Projection over index category "
690 << "requested, but projection data set doesn't contain index category" << endl ;
691 return frame ;
692 }
693
694 if (projectedVars.find(_indexCat.arg().GetName())) {
695 projIndex=true ;
696 }
697
698 } else {
699 // *** Error checking for a composite index category ***
700
701 // Determine if any servers of the index category are in the projectedVars
702 RooArgSet projIdxServers ;
703 bool anyServers(false) ;
704 for (const auto server : _indexCat->servers()) {
705 if (projectedVars.find(server->GetName())) {
706 anyServers=true ;
707 projIdxServers.add(*server) ;
708 }
709 }
710
711 // Check that the projection dataset contains all the
712 // index category components we're projecting over
713
714 // Determine if all projected servers of the index category are in the projection dataset
715 bool allServers(true) ;
716 std::string missing;
717 for (const auto server : projIdxServers) {
718 if (!projData->get()->find(server->GetName())) {
719 allServers=false ;
720 missing = server->GetName();
721 }
722 }
723
724 if (!allServers) {
725 coutE(Plotting) << "RooSimultaneous::plotOn(" << GetName()
726 << ") ERROR: Projection dataset doesn't contain complete set of index categories to do projection."
727 << "\n\tcategory " << missing << " is missing." << endl ;
728 return frame ;
729 }
730
731 if (anyServers) {
732 projIndex = true ;
733 }
734 }
735
736 // Calculate relative weight fractions of components
737 std::unique_ptr<Roo1DTable> wTable( projData->table(_indexCat.arg()) );
738
739 // Clone the index category to be able to cycle through the category states for plotting without
740 // affecting the category state of our instance
741 std::unique_ptr<RooArgSet> idxCloneSet( RooArgSet(*_indexCat).snapshot(true) );
742 auto idxCatClone = static_cast<RooAbsCategoryLValue*>( idxCloneSet->find(_indexCat->GetName()) );
743 assert(idxCatClone);
744
745 // Make list of category columns to exclude from projection data
746 std::unique_ptr<RooArgSet> idxCompSliceSet( idxCatClone->getObservables(frame->getNormVars()) );
747
748 // If we don't project over the index, just do the regular plotOn
749 if (!projIndex) {
750
751 coutI(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") plot on " << frame->getPlotVar()->GetName()
752 << " represents a slice in the index category (" << _indexCat.arg().GetName() << ")" << endl ;
753
754 // Reduce projData: take out fitCat (component) columns and entries that don't match selected slice
755 // Construct cut string to only select projection data event that match the current slice
756
757 // Make cut string to exclude rows from projection data
758 TString cutString ;
759 bool first(true) ;
760 for (const auto arg : *idxCompSliceSet) {
761 auto idxComp = static_cast<RooCategory*>(arg);
762 RooAbsArg* slicedComponent = nullptr;
763 if (sliceSet && (slicedComponent = sliceSet->find(*idxComp)) != nullptr) {
764 auto theCat = static_cast<const RooAbsCategory*>(slicedComponent);
765 idxComp->setIndex(theCat->getCurrentIndex(), false);
766 }
767
768 if (!first) {
769 cutString.Append("&&") ;
770 } else {
771 first=false ;
772 }
773 cutString.Append(Form("%s==%d",idxComp->GetName(),idxComp->getCurrentIndex())) ;
774 }
775
776 // Make temporary projData without RooSim index category components
777 RooArgSet projDataVars(*projData->get()) ;
778 projDataVars.remove(*idxCompSliceSet,true,true) ;
779
780 std::unique_ptr<RooAbsData> projDataTmp( const_cast<RooAbsData*>(projData)->reduce(projDataVars,cutString) );
781
782 // Override normalization and projection dataset
783 RooCmdArg tmp1 = RooFit::Normalization(scaleFactor*wTable->getFrac(idxCatClone->getCurrentLabel()),stype) ;
784 RooCmdArg tmp2 = RooFit::ProjWData(*projDataSet,*projDataTmp) ;
785
786 // WVE -- do not adjust normalization for asymmetry plots
787 RooLinkedList cmdList2(cmdList) ;
788 if (!cmdList.find("Asymmetry")) {
789 cmdList2.Add(&tmp1) ;
790 }
791 cmdList2.Add(&tmp2) ;
792
793 // Plot single component
794 RooPlot* retFrame = getPdf(idxCatClone->getCurrentLabel())->plotOn(frame,cmdList2);
795 return retFrame ;
796 }
797
798 // If we project over the index, plot using a temporary RooAddPdf
799 // using the weights from the data as coefficients
800
801 // Build the list of indexCat components that are sliced
802 idxCompSliceSet->remove(projectedVars,true,true) ;
803
804 // Make a new expression that is the weighted sum of requested components
805 RooArgList pdfCompList ;
806 RooArgList wgtCompList ;
807//RooAbsPdf* pdf ;
808 RooRealProxy* proxy ;
810 double sumWeight(0) ;
811 while((proxy=(RooRealProxy*)pIter.Next())) {
812
813 idxCatClone->setLabel(proxy->name()) ;
814
815 // Determine if this component is the current slice (if we slice)
816 bool skip(false) ;
817 for (const auto idxSliceCompArg : *idxCompSliceSet) {
818 const auto idxSliceComp = static_cast<RooAbsCategory*>(idxSliceCompArg);
819 RooAbsCategory* idxComp = (RooAbsCategory*) idxCloneSet->find(idxSliceComp->GetName()) ;
820 if (idxComp->getCurrentIndex()!=idxSliceComp->getCurrentIndex()) {
821 skip=true ;
822 break ;
823 }
824 }
825 if (skip) continue ;
826
827 // Instantiate a RRV holding this pdfs weight fraction
828 RooRealVar *wgtVar = new RooRealVar(proxy->name(),"coef",wTable->getFrac(proxy->name())) ;
829 wgtCompList.addOwned(*wgtVar) ;
830 sumWeight += wTable->getFrac(proxy->name()) ;
831
832 // Add the PDF to list list
833 pdfCompList.add(proxy->arg()) ;
834 }
835
836 TString plotVarName(GetName()) ;
837 RooAddPdf *plotVar = new RooAddPdf(plotVarName,"weighted sum of RS components",pdfCompList,wgtCompList) ;
838
839 // Fix appropriate coefficient normalization in plot function
840 if (_plotCoefNormSet.getSize()>0) {
842 }
843
844 std::unique_ptr<RooAbsData> projDataTmp;
845 RooArgSet projSetTmp ;
846 if (projData) {
847
848 // Construct cut string to only select projection data event that match the current slice
849 TString cutString ;
850 if (idxCompSliceSet->getSize()>0) {
851 bool first(true) ;
852 for (const auto idxSliceCompArg : *idxCompSliceSet) {
853 const auto idxSliceComp = static_cast<RooAbsCategory*>(idxSliceCompArg);
854 if (!first) {
855 cutString.Append("&&") ;
856 } else {
857 first=false ;
858 }
859 cutString.Append(Form("%s==%d",idxSliceComp->GetName(),idxSliceComp->getCurrentIndex())) ;
860 }
861 }
862
863 // Make temporary projData without RooSim index category components
864 RooArgSet projDataVars(*projData->get()) ;
865 RooArgSet* idxCatServers = _indexCat.arg().getObservables(frame->getNormVars()) ;
866
867 projDataVars.remove(*idxCatServers,true,true) ;
868
869 if (idxCompSliceSet->getSize()>0) {
870 projDataTmp.reset( const_cast<RooAbsData*>(projData)->reduce(projDataVars,cutString) );
871 } else {
872 projDataTmp.reset( const_cast<RooAbsData*>(projData)->reduce(projDataVars) );
873 }
874
875
876
877 if (projSet) {
878 projSetTmp.add(*projSet) ;
879 projSetTmp.remove(*idxCatServers,true,true);
880 }
881
882
883 delete idxCatServers ;
884 }
885
886
887 if (_indexCat.arg().isDerived() && idxCompSliceSet->getSize()>0) {
888 coutI(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") plot on " << frame->getPlotVar()->GetName()
889 << " represents a slice in index category components " << *idxCompSliceSet << endl ;
890
891 RooArgSet* idxCompProjSet = _indexCat.arg().getObservables(frame->getNormVars()) ;
892 idxCompProjSet->remove(*idxCompSliceSet,true,true) ;
893 if (idxCompProjSet->getSize()>0) {
894 coutI(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") plot on " << frame->getPlotVar()->GetName()
895 << " averages with data index category components " << *idxCompProjSet << endl ;
896 }
897 delete idxCompProjSet ;
898 } else {
899 coutI(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") plot on " << frame->getPlotVar()->GetName()
900 << " averages with data index category (" << _indexCat.arg().GetName() << ")" << endl ;
901 }
902
903
904 // Override normalization and projection dataset
905 RooLinkedList cmdList2(cmdList) ;
906
907 RooCmdArg tmp1 = RooFit::Normalization(scaleFactor*sumWeight,stype) ;
908 RooCmdArg tmp2 = RooFit::ProjWData(*projDataSet,*projDataTmp) ;
909 // WVE -- do not adjust normalization for asymmetry plots
910 if (!cmdList.find("Asymmetry")) {
911 cmdList2.Add(&tmp1) ;
912 }
913 cmdList2.Add(&tmp2) ;
914
915 RooPlot* frame2 ;
916 if (projSetTmp.getSize()>0) {
917 // Plot temporary function
918 RooCmdArg tmp3 = RooFit::Project(projSetTmp) ;
919 cmdList2.Add(&tmp3) ;
920 frame2 = plotVar->plotOn(frame,cmdList2) ;
921 } else {
922 // Plot temporary function
923 frame2 = plotVar->plotOn(frame,cmdList2) ;
924 }
925
926 // Cleanup
927 delete plotVar ;
928
929 return frame2 ;
930}
931
932
933
934////////////////////////////////////////////////////////////////////////////////
935/// OBSOLETE -- Retained for backward compatibility
936
937RooPlot* RooSimultaneous::plotOn(RooPlot *frame, Option_t* drawOptions, double scaleFactor,
938 ScaleType stype, const RooAbsData* projData, const RooArgSet* projSet,
939 double /*precision*/, bool /*shiftToZero*/, const RooArgSet* /*projDataSet*/,
940 double /*rangeLo*/, double /*rangeHi*/, RooCurve::WingMode /*wmode*/) const
941{
942 // Make command list
943 RooLinkedList cmdList ;
944 cmdList.Add(new RooCmdArg(RooFit::DrawOption(drawOptions))) ;
945 cmdList.Add(new RooCmdArg(RooFit::Normalization(scaleFactor,stype))) ;
946 if (projData) cmdList.Add(new RooCmdArg(RooFit::ProjWData(*projData))) ;
947 if (projSet) cmdList.Add(new RooCmdArg(RooFit::Project(*projSet))) ;
948
949 // Call new method
950 RooPlot* ret = plotOn(frame,cmdList) ;
951
952 // Cleanup
953 cmdList.Delete() ;
954 return ret ;
955}
956
957
958
959////////////////////////////////////////////////////////////////////////////////
960/// Interface function used by test statistics to freeze choice of observables
961/// for interpretation of fraction coefficients. Needed here because a RooSimultaneous
962/// works like a RooAddPdf when plotted
963
964void RooSimultaneous::selectNormalization(const RooArgSet* normSet, bool /*force*/)
965{
967 if (normSet) _plotCoefNormSet.add(*normSet) ;
968}
969
970
971////////////////////////////////////////////////////////////////////////////////
972/// Interface function used by test statistics to freeze choice of range
973/// for interpretation of fraction coefficients. Needed here because a RooSimultaneous
974/// works like a RooAddPdf when plotted
975
976void RooSimultaneous::selectNormalizationRange(const char* normRange2, bool /*force*/)
977{
978 _plotCoefNormRange = RooNameReg::ptr(normRange2) ;
979}
980
981
982
983
984////////////////////////////////////////////////////////////////////////////////
985
987 const RooArgSet* auxProto, bool verbose, bool autoBinned, const char* binnedTag) const
988{
989 const char* idxCatName = _indexCat.arg().GetName() ;
990
991 if (vars.find(idxCatName) && prototype==0
992 && (auxProto==0 || auxProto->empty())
993 && (autoBinned || (binnedTag && strlen(binnedTag)))) {
994
995 // Return special generator config that can also do binned generation for selected states
996 return new RooSimSplitGenContext(*this,vars,verbose,autoBinned,binnedTag) ;
997
998 } else {
999
1000 // Return regular generator config ;
1001 return genContext(vars,prototype,auxProto,verbose) ;
1002 }
1003}
1004
1005
1006
1007////////////////////////////////////////////////////////////////////////////////
1008/// Return specialized generator context for simultaneous p.d.f.s
1009
1011 const RooArgSet* auxProto, bool verbose) const
1012{
1013 const char* idxCatName = _indexCat.arg().GetName() ;
1014 const RooArgSet* protoVars = prototype ? prototype->get() : 0 ;
1015
1016 if (vars.find(idxCatName) || (protoVars && protoVars->find(idxCatName))) {
1017
1018 // Generating index category: return special sim-context
1019 return new RooSimGenContext(*this,vars,prototype,auxProto,verbose) ;
1020
1021 } else if (_indexCat.arg().isDerived()) {
1022 // Generating dependents of a derived index category
1023
1024 // Determine if we none,any or all servers
1025 bool anyServer(false), allServers(true) ;
1026 if (prototype) {
1027 TIterator* sIter = _indexCat.arg().serverIterator() ;
1028 RooAbsArg* server ;
1029 while((server=(RooAbsArg*)sIter->Next())) {
1030 if (prototype->get()->find(server->GetName())) {
1031 anyServer=true ;
1032 } else {
1033 allServers=false ;
1034 }
1035 }
1036 delete sIter ;
1037 } else {
1038 allServers=true ;
1039 }
1040
1041 if (allServers) {
1042 // Use simcontext if we have all servers
1043
1044 return new RooSimGenContext(*this,vars,prototype,auxProto,verbose) ;
1045 } else if (!allServers && anyServer) {
1046 // Abort if we have only part of the servers
1047 coutE(Plotting) << "RooSimultaneous::genContext: ERROR: prototype must include either all "
1048 << " components of the RooSimultaneous index category or none " << endl ;
1049 return 0 ;
1050 }
1051 // Otherwise make single gencontext for current state
1052 }
1053
1054 // Not generating index cat: return context for pdf associated with present index state
1056 if (!proxy) {
1057 coutE(InputArguments) << "RooSimultaneous::genContext(" << GetName()
1058 << ") ERROR: no PDF associated with current state ("
1059 << _indexCat.arg().GetName() << "=" << _indexCat.arg().getCurrentLabel() << ")" << endl ;
1060 return 0 ;
1061 }
1062 return ((RooAbsPdf*)proxy->absArg())->genContext(vars,prototype,auxProto,verbose) ;
1063}
1064
1065
1066
1067
1068////////////////////////////////////////////////////////////////////////////////
1069
1071 const RooArgSet* nset,
1072 double scaleFactor,
1073 bool correctForBinVolume,
1074 bool showProgress) const
1075{
1076 if (RooAbsReal::fillDataHist (hist, nset, scaleFactor,
1077 correctForBinVolume, showProgress) == 0)
1078 return 0;
1079
1080 const double sum = hist->sumEntries();
1081 if (sum != 0) {
1082 for (int i=0 ; i<hist->numEntries() ; i++) {
1083 hist->set(i, hist->weight(i) / sum, 0.);
1084 }
1085 }
1086
1087 return hist;
1088}
1089
1090
1091
1092
1093////////////////////////////////////////////////////////////////////////////////
1094/// Special generator interface for generation of 'global observables' -- for RooStats tools
1095
1097{
1098 // Make set with clone of variables (placeholder for output)
1099 RooArgSet* globClone = (RooArgSet*) whatVars.snapshot() ;
1100
1101 RooDataSet* data = new RooDataSet("gensimglobal","gensimglobal",whatVars) ;
1102
1103 for (Int_t i=0 ; i<nEvents ; i++) {
1104 for (const auto& nameIdx : indexCat()) {
1105
1106 // Get pdf associated with state from simpdf
1107 RooAbsPdf* pdftmp = getPdf(nameIdx.first.c_str());
1108
1109 // Generate only global variables defined by the pdf associated with this state
1110 RooArgSet* globtmp = pdftmp->getObservables(whatVars) ;
1111 RooDataSet* tmp = pdftmp->generate(*globtmp,1) ;
1112
1113 // Transfer values to output placeholder
1114 globClone->assign(*tmp->get(0)) ;
1115
1116 // Cleanup
1117 delete globtmp ;
1118 delete tmp ;
1119 }
1120 data->add(*globClone) ;
1121 }
1122
1123 delete globClone ;
1124 return data ;
1125}
1126
1127
1128/// Wraps the components of this RooSimultaneous in RooBinSamplingPdfs.
1129/// \param[in] data The dataset to be used in the eventual fit, used to figure
1130/// out the observables and whether the dataset is binned.
1131/// \param[in] precision Precision argument for all created RooBinSamplingPdfs.
1133
1134 if (precision < 0.) return;
1135
1136 RooArgSet newSamplingPdfs;
1137
1138 for (auto const &item : this->indexCat()) {
1139
1140 auto const &catName = item.first;
1141 auto &pdf = *this->getPdf(catName.c_str());
1142
1143 if (auto newSamplingPdf = RooBinSamplingPdf::create(pdf, data, precision)) {
1144 // Set the "ORIGNAME" attribute the indicate to
1145 // RooAbsArg::redirectServers() wich pdf should be replaced by this
1146 // RooBinSamplingPdf in the RooSimultaneous.
1147 newSamplingPdf->setAttribute(
1148 (std::string("ORIGNAME:") + pdf.GetName()).c_str());
1149 newSamplingPdfs.addOwned(std::move(newSamplingPdf));
1150 }
1151 }
1152
1153 this->redirectServers(newSamplingPdfs, false, true);
1154 this->addOwnedComponents(std::move(newSamplingPdfs));
1155}
1156
1157
1158/// Wraps the components of this RooSimultaneous in RooBinSamplingPdfs, with a
1159/// different precision parameter for each component.
1160/// \param[in] data The dataset to be used in the eventual fit, used to figure
1161/// out the observables and whether the dataset is binned.
1162/// \param[in] precisions The map that gives the precision argument for each
1163/// component in the RooSimultaneous. The keys are the pdf names. If
1164/// there is no value for a given component, it will not use the bin
1165/// integration. Otherwise, the value has the same meaning than in
1166/// the IntegrateBins() command argument for RooAbsPdf::fitTo().
1167/// \param[in] useCategoryNames If this flag is set, the category names will be
1168/// used to look up the precision in the precisions map instead of
1169/// the pdf names.
1171 std::map<std::string, double> const& precisions,
1172 bool useCategoryNames /*=false*/) {
1173
1174 constexpr double defaultPrecision = -1.;
1175
1176 RooArgSet newSamplingPdfs;
1177
1178 for (auto const &item : this->indexCat()) {
1179
1180 auto const &catName = item.first;
1181 auto &pdf = *this->getPdf(catName.c_str());
1182 std::string pdfName = pdf.GetName();
1183
1184 auto found = precisions.find(useCategoryNames ? catName : pdfName);
1185 const double precision =
1186 found != precisions.end() ? found->second : defaultPrecision;
1187 if (precision < 0.)
1188 continue;
1189
1190 if (auto newSamplingPdf = RooBinSamplingPdf::create(pdf, data, precision)) {
1191 // Set the "ORIGNAME" attribute the indicate to
1192 // RooAbsArg::redirectServers() wich pdf should be replaced by this
1193 // RooBinSamplingPdf in the RooSimultaneous.
1194 newSamplingPdf->setAttribute(
1195 (std::string("ORIGNAME:") + pdf.GetName()).c_str());
1196 newSamplingPdfs.addOwned(std::move(newSamplingPdf));
1197 }
1198 }
1199
1200 this->redirectServers(newSamplingPdfs, false, true);
1201 this->addOwnedComponents(std::move(newSamplingPdfs));
1202}
#define coutI(a)
Definition: RooMsgService.h:34
#define cxcoutD(a)
Definition: RooMsgService.h:85
#define coutW(a)
Definition: RooMsgService.h:36
#define coutE(a)
Definition: RooMsgService.h:37
RooTemplateProxy< RooAbsReal > RooRealProxy
Compatibility typedef replacing the old RooRealProxy class.
Definition: RooRealProxy.h:23
const char Option_t
Definition: RtypesCore.h:66
#define ClassImp(name)
Definition: Rtypes.h:375
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t type
char name[80]
Definition: TGX11.cxx:110
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition: TString.cxx:2447
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:72
friend class RooDataSet
Definition: RooAbsArg.h:690
RooArgSet * getObservables(const RooArgSet &set, bool valueOnly=true) const
Given a set of possible observables, return the observables that this PDF depends on.
Definition: RooAbsArg.h:312
friend class RooArgSet
Definition: RooAbsArg.h:645
bool addOwnedComponents(const RooAbsCollection &comps)
Take ownership of the contents of 'comps'.
Definition: RooAbsArg.cxx:2185
bool redirectServers(const RooAbsCollection &newServerList, bool mustReplaceAll=false, bool nameChange=false, bool isRecursionStep=false)
Replace all direct servers of this object with the new servers in newServerList.
Definition: RooAbsArg.cxx:1025
bool dependsOn(const RooAbsCollection &serverList, const RooAbsArg *ignoreArg=0, bool valueOnly=false) const
Test whether we depend on (ie, are served by) any object in the specified collection.
Definition: RooAbsArg.cxx:829
const RefCountList_t & servers() const
List of all servers of this object.
Definition: RooAbsArg.h:199
RooArgSet * getVariables(bool stripDisconnected=true) const
Return RooArgSet with all variables (tree leaf nodes of expresssion tree)
Definition: RooAbsArg.cxx:2057
virtual bool isDerived() const
Does value or shape of this arg depend on any other arg?
Definition: RooAbsArg.h:92
virtual bool isFundamental() const
Is this object a fundamental type that can be added to a dataset? Fundamental-type subclasses overrid...
Definition: RooAbsArg.h:242
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 printError=true)=0
Change category state by specifying a state name.
A space to attach TBranches.
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.
virtual bool addOwned(RooAbsArg &var, bool silent=false)
Add an argument and transfer the ownership to the collection.
virtual bool remove(const RooAbsArg &var, bool silent=false, bool matchByNameOnly=false)
Remove the specified argument from our list.
bool contains(const RooAbsArg &var) const
Check if collection contains an argument with the same name as var.
bool empty() const
Int_t getSize() const
Return the number of elements in the collection.
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
void assign(const RooAbsCollection &other) const
Sets the value, cache and constant attribute of any argument in our set that also appears in the othe...
Storage_t::size_type size() const
RooAbsArg * first() const
TObject * Clone(const char *newname=0) const override
Make a clone of an object using the Streamer facility.
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:62
virtual const RooArgSet * get() const
Definition: RooAbsData.h:106
virtual Roo1DTable * table(const RooArgSet &catSet, const char *cuts="", const char *opts="") const
Construct table for product of categories in catSet.
Definition: RooAbsData.cxx:889
RooAbsGenContext is the abstract base class for generator contexts of RooAbsPdf objects.
friend class CacheElem
The cache manager.
Definition: RooAbsPdf.h:367
bool mustBeExtended() const
If true PDF must provide extended likelihood term.
Definition: RooAbsPdf.h:267
RooArgSet const * _normSet
Normalization integral (owned by _normMgr)
Definition: RooAbsPdf.h:355
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:60
bool canBeExtended() const
If true, PDF can provide extended likelihood term.
Definition: RooAbsPdf.h:263
@ CanBeExtended
Definition: RooAbsPdf.h:257
@ MustBeExtended
Definition: RooAbsPdf.h:257
@ CanNotBeExtended
Definition: RooAbsPdf.h:257
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 override
Helper calling plotOn(RooPlot*, RooLinkedList&) const.
Definition: RooAbsPdf.h:126
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:64
RooDataHist * fillDataHist(RooDataHist *hist, const RooArgSet *nset, double scaleFactor, bool correctForBinVolume=false, bool showProgress=false) const
Fill a RooDataHist with values sampled from this function at the bin centers.
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:94
virtual void fixAddCoefNormalization(const RooArgSet &addNormSet=RooArgSet(), bool force=true)
Fix the interpretation of the coefficient of any RooAddPdf component in the expression tree headed by...
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:553
bool plotSanityChecks(RooPlot *frame) const
Utility function for plotOn(), perform general sanity check on frame to ensure safe plotting operatio...
friend class RooAddPdf
Definition: RooAbsReal.h:558
void makeProjectionSet(const RooAbsArg *plotVar, const RooArgSet *allVars, RooArgSet &projectedVars, bool silent) const
Utility function for plotOn() that constructs the set of observables to project when plotting ourselv...
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:57
RooArgSet * snapshot(bool deepCopy=true) const
Use RooAbsCollection::snapshot(), but return as RooArgSet.
Definition: RooArgSet.h:180
static std::unique_ptr< RooAbsPdf > create(RooAbsPdf &pdf, RooAbsData const &data, double precision)
Creates a wrapping RooBinSamplingPdf if appropriate.
T * getObj(const RooArgSet *nset, Int_t *sterileIndex=0, const TNamed *isetRangeName=0)
Getter function without integration set.
T * getObjByIndex(Int_t index) const
Retrieve payload object by slot index.
Int_t lastIndex() const
Return index of slot used in last get or set operation.
Int_t setObj(const RooArgSet *nset, T *obj, const TNamed *isetRangeName=0)
Setter function without integration set.
RooCategory is an object to represent discrete states.
Definition: RooCategory.h:28
RooCmdArg is a named container for two doubles, two integers two object points and three string point...
Definition: RooCmdArg.h:26
Class RooCmdConfig is a configurable parser for RooCmdArg named arguments.
Definition: RooCmdConfig.h:31
void removeAll() override
Remove all argument inset using remove(const RooAbsArg&).
bool add(const RooAbsArg &var, bool valueServer, bool shapeServer, bool silent)
Overloaded RooCollection_t::add() method insert object into set and registers object as server to own...
The RooDataHist is a container class to hold N-dimensional binned data.
Definition: RooDataHist.h:45
double weight(std::size_t i) const
Return weight of i-th bin.
Definition: RooDataHist.h:110
void set(std::size_t binNumber, double weight, double wgtErr)
Set bin content of bin that was last loaded with get(std::size_t).
Int_t numEntries() const override
Return the number of bins.
double sumEntries() const override
Sum the weights of all bins.
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:55
const RooArgSet * get(Int_t index) const override
Return RooArgSet with coordinates of event 'index'.
RooLinkedList is an collection class for internal use, storing a collection of RooAbsArg pointers in ...
Definition: RooLinkedList.h:38
void Delete(Option_t *o=0) override
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:81
A RooPlot is a plot frame and a container for graphics objects within that frame.
Definition: RooPlot.h:43
const RooArgSet * getNormVars() const
Definition: RooPlot.h:149
void remove(const char *name=0, bool deleteToo=true)
Remove object with given name, or last object added if no name is given.
Definition: RooPlot.cxx:913
RooAbsRealLValue * getPlotVar() const
Definition: RooPlot.h:140
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:40
RooSimultaneous facilitates simultaneous fitting of multiple PDFs to subsets of a given dataset.
double evaluate() const override
Return the current value: the value of the PDF associated with the current index category state.
Int_t _numPdf
Number of registered PDFs.
TList _pdfProxyList
List of PDF proxies (named after applicable category state)
RooObjCacheManager _partIntMgr
! Component normalization manager
void selectNormalization(const RooArgSet *depSet=0, bool force=false) override
Interface function used by test statistics to freeze choice of observables for interpretation of frac...
RooAbsGenContext * autoGenContext(const RooArgSet &vars, const RooDataSet *prototype=0, const RooArgSet *auxProto=0, bool verbose=false, bool autoBinned=true, const char *binnedTag="") const override
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 override
RooDataSet * generateSimGlobal(const RooArgSet &whatVars, Int_t nEvents) override
Special generator interface for generation of 'global observables' – for RooStats tools.
Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &numVars, const RooArgSet *normSet, const char *rangeName=0) const override
Forward determination of analytical integration capabilities to component p.d.f.s A unique code is as...
~RooSimultaneous() override
Destructor.
ExtendMode extendMode() const override
Examine the pdf components and check if one of them can be extended or must be extended It is enough ...
RooAbsGenContext * genContext(const RooArgSet &vars, const RooDataSet *prototype=0, const RooArgSet *auxProto=0, bool verbose=false) const override
Return specialized generator context for simultaneous p.d.f.s.
RooCategoryProxy _indexCat
Index category.
void selectNormalizationRange(const char *rangeName=0, bool force=false) override
Interface function used by test statistics to freeze choice of range for interpretation of fraction c...
double expectedEvents(const RooArgSet *nset) const override
Return the number of expected events: If the index is in nset, then return the sum of the expected ev...
friend class RooSimGenContext
void initialize(RooAbsCategoryLValue &inIndexCat, std::map< std::string, RooAbsPdf * > pdfMap)
double analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=0) const override
Return analytical integration defined by given code.
void wrapPdfsInBinSamplingPdfs(RooAbsData const &data, double precision)
Wraps the components of this RooSimultaneous in RooBinSamplingPdfs.
const TNamed * _plotCoefNormRange
bool addPdf(const RooAbsPdf &pdf, const char *catLabel)
Associate given PDF with index category state label 'catLabel'.
RooSetProxy _plotCoefNormSet
const RooAbsCategoryLValue & indexCat() const
friend class RooSimSplitGenContext
virtual RooDataHist * fillDataHist(RooDataHist *hist, const RooArgSet *nset, double scaleFactor, bool correctForBinVolume=false, bool showProgress=false) const
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.
bool setLabel(const char *label, bool printError=true) 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
TObject * FindObject(const char *name) const override
Find an object in this list using its name.
Definition: TList.cxx:578
void Add(TObject *obj) override
Definition: TList.h:81
TIterator * MakeIterator(Bool_t dir=kIterForward) const override
Return a list iterator.
Definition: TList.cxx:722
void Delete(Option_t *option="") override
Remove all objects from the list AND delete all heap based objects.
Definition: TList.cxx:470
TObject * At(Int_t idx) const override
Returns the object at position idx. Returns 0 if idx is out of range.
Definition: TList.cxx:357
const char * GetName() const override
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 binData=false)
RooCmdArg DrawOption(const char *opt)
RooCmdArg Project(const RooArgSet &projSet)
RooCmdArg Normalization(double 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:64
static constexpr double pc
Definition: first.py:1
const RooAbsCategoryLValue * subIndex
static uint64_t sum(uint64_t i)
Definition: Factory.cxx:2345