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,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 "RooRealVar.h"
53#include "RooAddPdf.h"
54#include "RooAbsData.h"
55#include "Roo1DTable.h"
56#include "RooSimGenContext.h"
58#include "RooDataSet.h"
59#include "RooCmdConfig.h"
60#include "RooNameReg.h"
61#include "RooGlobalFunc.h"
62#include "RooMsgService.h"
63#include "RooCategory.h"
64#include "RooSuperCategory.h"
65#include "RooDataHist.h"
66#include "RooRandom.h"
67#include "RooArgSet.h"
68#include "RooBinSamplingPdf.h"
69
70#include "ROOT/StringUtils.hxx"
71
72#include <iostream>
73
74namespace {
75
76std::map<std::string, RooAbsPdf *> createPdfMap(const RooArgList &inPdfList, RooAbsCategoryLValue &inIndexCat)
77{
78 std::map<std::string, RooAbsPdf *> pdfMap;
79 auto indexCatIt = inIndexCat.begin();
80 for (unsigned int i = 0; i < inPdfList.size(); ++i) {
81 auto pdf = static_cast<RooAbsPdf *>(&inPdfList[i]);
82 const auto &nameIdx = (*indexCatIt++);
83 pdfMap[nameIdx.first] = pdf;
84 }
85 return pdfMap;
86}
87
88} // namespace
89
91
92void RooSimultaneous::InitializationOutput::addPdf(const RooAbsPdf &pdf, std::string const &catLabel)
93{
94 finalPdfs.push_back(&pdf);
95 finalCatLabels.emplace_back(catLabel);
96}
97
98using namespace std;
99
101
102
103
104////////////////////////////////////////////////////////////////////////////////
105/// Constructor with index category. PDFs associated with indexCat
106/// states can be added after construction with the addPdf() function.
107///
108/// RooSimultaneous can function without having a PDF associated
109/// with every single state. The normalization in such cases is taken
110/// from the number of registered PDFs, but getVal() will assert if
111/// when called for an unregistered index state.
112
113RooSimultaneous::RooSimultaneous(const char *name, const char *title,
114 RooAbsCategoryLValue& inIndexCat) :
115 RooSimultaneous{name, title, std::map<std::string, RooAbsPdf*>{}, inIndexCat}
116{
117}
118
119
120////////////////////////////////////////////////////////////////////////////////
121/// Constructor from index category and full list of PDFs.
122/// In this constructor form, a PDF must be supplied for each indexCat state
123/// to avoid ambiguities. The PDFs are associated with the states of the
124/// index category as they appear when iterating through the category states
125/// with RooAbsCategory::begin() and RooAbsCategory::end(). This usually means
126/// they are associated by ascending index numbers.
127///
128/// PDFs may not overlap (i.e. share any variables) with the index category (function)
129
130RooSimultaneous::RooSimultaneous(const char *name, const char *title,
131 const RooArgList& inPdfList, RooAbsCategoryLValue& inIndexCat) :
132 RooSimultaneous{name, title, createPdfMap(inPdfList, inIndexCat), inIndexCat}
133{
134 if (inPdfList.size() != inIndexCat.size()) {
135 std::stringstream errMsg;
136 errMsg << "RooSimultaneous::ctor(" << GetName()
137 << " ERROR: Number PDF list entries must match number of index category states, no PDFs added";
138 coutE(InputArguments) << errMsg.str() << std::endl;
139 throw std::invalid_argument(errMsg.str());
140 }
141}
142
143
144////////////////////////////////////////////////////////////////////////////////
145
146RooSimultaneous::RooSimultaneous(const char *name, const char *title, std::map<string, RooAbsPdf *> pdfMap,
147 RooAbsCategoryLValue &inIndexCat)
148 : RooSimultaneous(name, title, std::move(*initialize(name ? name : "", inIndexCat, pdfMap)))
149{
150}
151
153 : RooAbsPdf(name, title),
154 _plotCoefNormSet("!plotCoefNormSet", "plotCoefNormSet", this, false, false),
155 _partIntMgr(this, 10),
156 _indexCat("indexCat", "Index category", this, *initInfo.indexCat)
157{
158 for (std::size_t i = 0; i < initInfo.finalPdfs.size(); ++i) {
159 addPdf(*initInfo.finalPdfs[i], initInfo.finalCatLabels[i].c_str());
160 }
161
162 // Take ownership of eventual super category
163 if (initInfo.superIndex) {
164 addOwnedComponents(std::move(initInfo.superIndex));
165 }
166}
167
168
169// This class cannot be locally defined in initialize as it cannot be
170// used as a template argument in that case
172 struct CompInfo {
176 std::unique_ptr<RooArgSet> subIndexComps;
177 } ;
178}
179
180std::unique_ptr<RooSimultaneous::InitializationOutput>
182 std::map<std::string, RooAbsPdf *> const& pdfMap)
183
184{
185 auto out = std::make_unique<RooSimultaneous::InitializationOutput>();
186 out->indexCat = &inIndexCat;
187
188 // First see if there are any RooSimultaneous input components
189 bool simComps(false) ;
190 for (auto const& item : pdfMap) {
191 if (dynamic_cast<RooSimultaneous*>(item.second)) {
192 simComps = true ;
193 break ;
194 }
195 }
196
197 // If there are no simultaneous component p.d.f. do simple processing through addPdf()
198 if (!simComps) {
199 for (auto const& item : pdfMap) {
200 out->addPdf(*item.second,item.first);
201 }
202 return out;
203 }
204
205 std::string msgPrefix = "RooSimultaneous::initialize(" + name + ") ";
206
207 // Issue info message that we are about to do some rearraning
208 oocoutI(nullptr, InputArguments) << msgPrefix << "INFO: one or more input component of simultaneous p.d.f.s are"
209 << " simultaneous p.d.f.s themselves, rewriting composite expressions as one-level simultaneous p.d.f. in terms of"
210 << " final constituents and extended index category" << std::endl;
211
212
213 RooArgSet allAuxCats ;
214 std::map<string,RooSimultaneousAux::CompInfo> compMap ;
215 for (auto const& item : pdfMap) {
217 ci.pdf = item.second ;
218 RooSimultaneous* simComp = dynamic_cast<RooSimultaneous*>(item.second) ;
219 if (simComp) {
220 ci.simPdf = simComp ;
221 ci.subIndex = &simComp->indexCat() ;
222 ci.subIndexComps = simComp->indexCat().isFundamental()
223 ? std::make_unique<RooArgSet>(simComp->indexCat())
224 : std::unique_ptr<RooArgSet>(simComp->indexCat().getVariables());
225 allAuxCats.add(*ci.subIndexComps,true) ;
226 } else {
227 ci.simPdf = nullptr;
228 ci.subIndex = nullptr;
229 }
230 compMap[item.first] = std::move(ci);
231 }
232
233 // Construct the 'superIndex' from the nominal index category and all auxiliary components
234 RooArgSet allCats(inIndexCat) ;
235 allCats.add(allAuxCats) ;
236 std::string siname = name + "_index";
237 out->superIndex = std::make_unique<RooSuperCategory>(siname.c_str(),siname.c_str(),allCats) ;
238 auto *superIndex = out->superIndex.get();
239 out->indexCat = superIndex;
240
241 // Now process each of original pdf/state map entries
242 for (auto const& citem : compMap) {
243
244 RooArgSet repliCats(allAuxCats) ;
245 if (citem.second.subIndexComps) {
246 repliCats.remove(*citem.second.subIndexComps) ;
247 }
248 inIndexCat.setLabel(citem.first.c_str()) ;
249
250 if (!citem.second.simPdf) {
251
252 // Entry is a plain p.d.f. assign it to every state permutation of the repliCats set
253 RooSuperCategory repliSuperCat("tmp","tmp",repliCats) ;
254
255 // Iterator over all states of repliSuperCat
256 for (const auto& nameIdx : repliSuperCat) {
257 // Set value
258 repliSuperCat.setLabel(nameIdx.first) ;
259 // Retrieve corresponding label of superIndex
260 string superLabel = superIndex->getCurrentLabel() ;
261 out->addPdf(*citem.second.pdf,superLabel);
262 oocxcoutD(static_cast<RooAbsArg*>(nullptr), InputArguments) << msgPrefix
263 << "assigning pdf " << citem.second.pdf->GetName() << " to super label " << superLabel << endl ;
264 }
265 } else {
266
267 // Entry is a simultaneous p.d.f
268
269 if (repliCats.empty()) {
270
271 // Case 1 -- No replication of components of RooSim component are required
272
273 for (const auto& type : *citem.second.subIndex) {
274 const_cast<RooAbsCategoryLValue*>(citem.second.subIndex)->setLabel(type.first.c_str());
275 string superLabel = superIndex->getCurrentLabel() ;
276 RooAbsPdf* compPdf = citem.second.simPdf->getPdf(type.first);
277 if (compPdf) {
278 out->addPdf(*compPdf,superLabel);
279 oocxcoutD(static_cast<RooAbsArg*>(nullptr), InputArguments) << msgPrefix
280 << "assigning pdf " << compPdf->GetName() << "(member of " << citem.second.pdf->GetName()
281 << ") to super label " << superLabel << endl ;
282 } else {
283 oocoutW(nullptr, InputArguments) << msgPrefix << "WARNING: No p.d.f. associated with label "
284 << type.second << " for component RooSimultaneous p.d.f " << citem.second.pdf->GetName()
285 << "which is associated with master index label " << citem.first << endl ;
286 }
287 }
288
289 } else {
290
291 // Case 2 -- Replication of components of RooSim component are required
292
293 // Make replication supercat
294 RooSuperCategory repliSuperCat("tmp","tmp",repliCats) ;
295
296 for (const auto& stype : *citem.second.subIndex) {
297 const_cast<RooAbsCategoryLValue*>(citem.second.subIndex)->setLabel(stype.first.c_str());
298
299 for (const auto& nameIdx : repliSuperCat) {
300 repliSuperCat.setLabel(nameIdx.first) ;
301 const string superLabel = superIndex->getCurrentLabel() ;
302 RooAbsPdf* compPdf = citem.second.simPdf->getPdf(stype.first);
303 if (compPdf) {
304 out->addPdf(*compPdf,superLabel);
305 oocxcoutD(static_cast<RooAbsArg*>(nullptr), InputArguments) << msgPrefix
306 << "assigning pdf " << compPdf->GetName() << "(member of " << citem.second.pdf->GetName()
307 << ") to super label " << superLabel << endl ;
308 } else {
309 oocoutW(nullptr, InputArguments) << msgPrefix << "WARNING: No p.d.f. associated with label "
310 << stype.second << " for component RooSimultaneous p.d.f " << citem.second.pdf->GetName()
311 << "which is associated with master index label " << citem.first << endl ;
312 }
313 }
314 }
315 }
316 }
317 }
318
319 return out;
320}
321
322
323////////////////////////////////////////////////////////////////////////////////
324/// Copy constructor
325
327 RooAbsPdf(other,name),
328 _plotCoefNormSet("!plotCoefNormSet",this,other._plotCoefNormSet),
329 _plotCoefNormRange(other._plotCoefNormRange),
330 _partIntMgr(other._partIntMgr,this),
331 _indexCat("indexCat",this,other._indexCat),
332 _numPdf(other._numPdf)
333{
334 // Copy proxy list
335 for(auto* proxy : static_range_cast<RooRealProxy*>(other._pdfProxyList)) {
336 _pdfProxyList.Add(new RooRealProxy(proxy->GetName(),this,*proxy)) ;
337 }
338}
339
340
341
342////////////////////////////////////////////////////////////////////////////////
343/// Destructor
344
346{
348}
349
350
351
352////////////////////////////////////////////////////////////////////////////////
353/// Return the p.d.f associated with the given index category name
354
356{
357 RooRealProxy* proxy = static_cast<RooRealProxy*>(_pdfProxyList.FindObject(catName));
358 return proxy ? static_cast<RooAbsPdf*>(proxy->absArg()) : nullptr;
359}
360
361
362
363////////////////////////////////////////////////////////////////////////////////
364/// Associate given PDF with index category state label 'catLabel'.
365/// The name state must be already defined in the index category.
366///
367/// RooSimultaneous can function without having a PDF associated
368/// with every single state. The normalization in such cases is taken
369/// from the number of registered PDFs, but getVal() will fail if
370/// called for an unregistered index state.
371///
372/// PDFs may not overlap (i.e. share any variables) with the index category (function).
373/// \param[in] pdf PDF to be added.
374/// \param[in] catLabel Name of the category state to be associated to the PDF.
375/// \return `true` in case of failure.
376
377bool RooSimultaneous::addPdf(const RooAbsPdf& pdf, const char* catLabel)
378{
379 // PDFs cannot overlap with the index category
380 if (pdf.dependsOn(_indexCat.arg())) {
381 coutE(InputArguments) << "RooSimultaneous::addPdf(" << GetName() << "): PDF '" << pdf.GetName()
382 << "' overlaps with index category '" << _indexCat.arg().GetName() << "'."<< endl ;
383 return true ;
384 }
385
386 // Each index state can only have one PDF associated with it
387 if (_pdfProxyList.FindObject(catLabel)) {
388 coutE(InputArguments) << "RooSimultaneous::addPdf(" << GetName() << "): index state '"
389 << catLabel << "' has already an associated PDF." << endl ;
390 return true ;
391 }
392
393 const RooSimultaneous* simPdf = dynamic_cast<const RooSimultaneous*>(&pdf) ;
394 if (simPdf) {
395
396 coutE(InputArguments) << "RooSimultaneous::addPdf(" << GetName()
397 << ") ERROR: you cannot add a RooSimultaneous component to a RooSimultaneous using addPdf()."
398 << " Use the constructor with RooArgList if input p.d.f.s or the map<string,RooAbsPdf&> instead." << endl ;
399 return true ;
400
401 } else {
402
403 // Create a proxy named after the associated index state
404 TObject* proxy = new RooRealProxy(catLabel,catLabel,this,const_cast<RooAbsPdf&>(pdf));
405 _pdfProxyList.Add(proxy) ;
406 _numPdf += 1 ;
407 }
408
409 return false ;
410}
411
412
413
414
415
416////////////////////////////////////////////////////////////////////////////////
417/// Examine the pdf components and check if one of them can be extended or must be extended
418/// It is enough to have one component that can be exteded or must be extended to return the flag in
419/// the total simultaneous pdf
420
422{
423 bool anyCanExtend(false) ;
424 bool anyMustExtend(false) ;
425
426 for (Int_t i=0 ; i<_numPdf ; i++) {
428 if (proxy) {
429 RooAbsPdf* pdf = (RooAbsPdf*) proxy->absArg() ;
430 //cout << " now processing pdf " << pdf->GetName() << endl;
431 if (pdf->canBeExtended()) {
432 //cout << "RooSim::extendedMode(" << GetName() << ") component " << pdf->GetName() << " can be extended"
433 // << endl;
434 anyCanExtend = true;
435 }
436 if (pdf->mustBeExtended()) {
437 //cout << "RooSim::extendedMode(" << GetName() << ") component " << pdf->GetName() << " MUST be extended" << endl;
438 anyMustExtend = true;
439 }
440 }
441 }
442 if (anyMustExtend) {
443 //cout << "RooSim::extendedMode(" << GetName() << ") returning MustBeExtended" << endl ;
444 return MustBeExtended ;
445 }
446 if (anyCanExtend) {
447 //cout << "RooSim::extendedMode(" << GetName() << ") returning CanBeExtended" << endl ;
448 return CanBeExtended ;
449 }
450 //cout << "RooSim::extendedMode(" << GetName() << ") returning CanNotBeExtended" << endl ;
451 return CanNotBeExtended ;
452}
453
454
455
456
457////////////////////////////////////////////////////////////////////////////////
458/// Return the current value:
459/// the value of the PDF associated with the current index category state
460
462{
463 // Retrieve the proxy by index name
465
466 //assert(proxy!=0) ;
467 if (proxy==0) return 0 ;
468
469 // Calculate relative weighting factor for sim-pdfs of all extendable components
470 double catFrac(1) ;
471 if (canBeExtended()) {
472 double nEvtCat = ((RooAbsPdf*)(proxy->absArg()))->expectedEvents(_normSet) ;
473
474 double nEvtTot(0) ;
475 for(auto * proxy2 : static_range_cast<RooRealProxy*>(_pdfProxyList)) {
476 nEvtTot += ((RooAbsPdf*)(proxy2->absArg()))->expectedEvents(_normSet) ;
477 }
478 catFrac=nEvtCat/nEvtTot ;
479 }
480
481 // Return the selected PDF value, normalized by the number of index states
482 return ((RooAbsPdf*)(proxy->absArg()))->getVal(_normSet)*catFrac ;
483}
484
485
486
487////////////////////////////////////////////////////////////////////////////////
488/// Return the number of expected events: If the index is in nset,
489/// then return the sum of the expected events of all components,
490/// otherwise return the number of expected events of the PDF
491/// associated with the current index category state
492
494{
495 if (nset->contains(_indexCat.arg())) {
496
497 double sum(0) ;
498
499 for(auto * proxy : static_range_cast<RooRealProxy*>(_pdfProxyList)) {
500 sum += ((RooAbsPdf*)(proxy->absArg()))->expectedEvents(nset) ;
501 }
502
503 return sum ;
504
505 } else {
506
507 // Retrieve the proxy by index name
509
510 //assert(proxy!=0) ;
511 if (proxy==0) return 0 ;
512
513 // Return the selected PDF value, normalized by the number of index states
514 return ((RooAbsPdf*)(proxy->absArg()))->expectedEvents(nset);
515 }
516}
517
518
519
520////////////////////////////////////////////////////////////////////////////////
521/// Forward determination of analytical integration capabilities to component p.d.f.s
522/// A unique code is assigned to the combined integration capabilities of all associated
523/// p.d.f.s
524
526 const RooArgSet* normSet, const char* rangeName) const
527{
528 // Declare that we can analytically integrate all requested observables
529 analVars.add(allVars) ;
530
531 // Retrieve (or create) the required partial integral list
532 Int_t code ;
533
534 // Check if this configuration was created before
535 CacheElem* cache = (CacheElem*) _partIntMgr.getObj(normSet,&analVars,0,RooNameReg::ptr(rangeName)) ;
536 if (cache) {
537 code = _partIntMgr.lastIndex() ;
538 return code+1 ;
539 }
540 cache = new CacheElem ;
541
542 // Create the partial integral set for this request
543 for(auto * proxy : static_range_cast<RooRealProxy*>(_pdfProxyList)) {
544 cache->_partIntList.addOwned(std::unique_ptr<RooAbsReal>{proxy->arg().createIntegral(analVars,normSet,0,rangeName)});
545 }
546
547 // Store the partial integral list and return the assigned code ;
548 code = _partIntMgr.setObj(normSet,&analVars,cache,RooNameReg::ptr(rangeName)) ;
549
550 return code+1 ;
551}
552
553
554
555////////////////////////////////////////////////////////////////////////////////
556/// Return analytical integration defined by given code
557
558double RooSimultaneous::analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* /*rangeName*/) const
559{
560 // No integration scenario
561 if (code==0) {
562 return getVal(normSet) ;
563 }
564
565 // Partial integration scenarios, rangeName already encoded in 'code'
566 CacheElem* cache = (CacheElem*) _partIntMgr.getObjByIndex(code-1) ;
567
569 Int_t idx = _pdfProxyList.IndexOf(proxy) ;
570 return ((RooAbsReal*)cache->_partIntList.at(idx))->getVal(normSet) ;
571}
572
573
574
575
576
577
578////////////////////////////////////////////////////////////////////////////////
579/// Back-end for plotOn() implementation on RooSimultaneous which
580/// needs special handling because a RooSimultaneous PDF cannot
581/// project out its index category via integration. plotOn() will
582/// abort if this is requested without providing a projection dataset.
583
585{
586 // Sanity checks
587 if (plotSanityChecks(frame)) return frame ;
588
589 // Extract projection configuration from command list
590 RooCmdConfig pc(Form("RooSimultaneous::plotOn(%s)",GetName())) ;
591 pc.defineString("sliceCatState","SliceCat",0,"",true) ;
592 pc.defineDouble("scaleFactor","Normalization",0,1.0) ;
593 pc.defineInt("scaleType","Normalization",0,RooAbsPdf::Relative) ;
594 pc.defineObject("sliceCatList","SliceCat",0,0,true) ;
595 // This dummy is needed for plotOn to recognize the "SliceCatMany" command.
596 // It is not used directly, but the "SliceCat" commands are nested in it.
597 // Removing this dummy definition results in "ERROR: unrecognized command: SliceCatMany".
598 pc.defineObject("dummy1","SliceCatMany",0) ;
599 pc.defineSet("projSet","Project",0) ;
600 pc.defineSet("sliceSet","SliceVars",0) ;
601 pc.defineSet("projDataSet","ProjData",0) ;
602 pc.defineObject("projData","ProjData",1) ;
603 pc.defineMutex("Project","SliceVars") ;
604 pc.allowUndefined() ; // there may be commands we don't handle here
605
606 // Process and check varargs
607 pc.process(cmdList) ;
608 if (!pc.ok(true)) {
609 return frame ;
610 }
611
612 const RooAbsData* projData = (const RooAbsData*) pc.getObject("projData") ;
613 const RooArgSet* projDataSet = pc.getSet("projDataSet");
614 const RooArgSet* sliceSetTmp = pc.getSet("sliceSet") ;
615 std::unique_ptr<RooArgSet> sliceSet( sliceSetTmp ? ((RooArgSet*) sliceSetTmp->Clone()) : nullptr );
616 const RooArgSet* projSet = pc.getSet("projSet") ;
617 double scaleFactor = pc.getDouble("scaleFactor") ;
618 ScaleType stype = (ScaleType) pc.getInt("scaleType") ;
619
620
621 // Look for category slice arguments and add them to the master slice list if found
622 const char* sliceCatState = pc.getString("sliceCatState",0,true) ;
623 const RooLinkedList& sliceCatList = pc.getObjectList("sliceCatList") ;
624 if (sliceCatState) {
625
626 // Make the master slice set if it doesnt exist
627 if (!sliceSet) {
628 sliceSet = std::make_unique<RooArgSet>();
629 }
630
631 // Prepare comma separated label list for parsing
632 auto catTokens = ROOT::Split(sliceCatState, ",");
633
634 // Loop over all categories provided by (multiple) Slice() arguments
635 unsigned int tokenIndex = 0;
636 for(auto * scat : static_range_cast<RooCategory*>(sliceCatList)) {
637 const char* slabel = tokenIndex >= catTokens.size() ? nullptr : catTokens[tokenIndex++].c_str();
638
639 if (slabel) {
640 // Set the slice position to the value indicated by slabel
641 scat->setLabel(slabel) ;
642 // Add the slice category to the master slice set
643 sliceSet->add(*scat,false) ;
644 }
645 }
646 }
647
648 // Check if we have a projection dataset
649 if (!projData) {
650 coutE(InputArguments) << "RooSimultaneous::plotOn(" << GetName() << ") ERROR: must have a projection dataset for index category" << endl ;
651 return frame ;
652 }
653
654 // Make list of variables to be projected
655 RooArgSet projectedVars ;
656 if (sliceSet) {
657 makeProjectionSet(frame->getPlotVar(),frame->getNormVars(),projectedVars,true) ;
658
659 // Take out the sliced variables
660 for (const auto sliceArg : *sliceSet) {
661 RooAbsArg* arg = projectedVars.find(sliceArg->GetName()) ;
662 if (arg) {
663 projectedVars.remove(*arg) ;
664 } else {
665 coutI(Plotting) << "RooAbsReal::plotOn(" << GetName() << ") slice variable "
666 << sliceArg->GetName() << " was not projected anyway" << endl ;
667 }
668 }
669 } else if (projSet) {
670 makeProjectionSet(frame->getPlotVar(),projSet,projectedVars,false) ;
671 } else {
672 makeProjectionSet(frame->getPlotVar(),frame->getNormVars(),projectedVars,true) ;
673 }
674
675 bool projIndex(false) ;
676
677 if (!_indexCat.arg().isDerived()) {
678 // *** Error checking for a fundamental index category ***
679 //cout << "RooSim::plotOn: index is fundamental" << endl ;
680
681 // Check that the provided projection dataset contains our index variable
682 if (!projData->get()->find(_indexCat.arg().GetName())) {
683 coutE(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") ERROR: Projection over index category "
684 << "requested, but projection data set doesn't contain index category" << endl ;
685 return frame ;
686 }
687
688 if (projectedVars.find(_indexCat.arg().GetName())) {
689 projIndex=true ;
690 }
691
692 } else {
693 // *** Error checking for a composite index category ***
694
695 // Determine if any servers of the index category are in the projectedVars
696 RooArgSet projIdxServers ;
697 bool anyServers(false) ;
698 for (const auto server : flattenedCatList()) {
699 if (projectedVars.find(server->GetName())) {
700 anyServers=true ;
701 projIdxServers.add(*server) ;
702 }
703 }
704
705 // Check that the projection dataset contains all the
706 // index category components we're projecting over
707
708 // Determine if all projected servers of the index category are in the projection dataset
709 bool allServers(true) ;
710 std::string missing;
711 for (const auto server : projIdxServers) {
712 if (!projData->get()->find(server->GetName())) {
713 allServers=false ;
714 missing = server->GetName();
715 }
716 }
717
718 if (!allServers) {
719 coutE(Plotting) << "RooSimultaneous::plotOn(" << GetName()
720 << ") ERROR: Projection dataset doesn't contain complete set of index categories to do projection."
721 << "\n\tcategory " << missing << " is missing." << endl ;
722 return frame ;
723 }
724
725 if (anyServers) {
726 projIndex = true ;
727 }
728 }
729
730 // Calculate relative weight fractions of components
731 std::unique_ptr<Roo1DTable> wTable( projData->table(_indexCat.arg()) );
732
733 // Clone the index category to be able to cycle through the category states for plotting without
734 // affecting the category state of our instance
735 RooArgSet idxCloneSet;
736 RooArgSet(*_indexCat).snapshot(idxCloneSet, true);
737 auto idxCatClone = static_cast<RooAbsCategoryLValue*>( idxCloneSet.find(_indexCat->GetName()) );
738 assert(idxCatClone);
739
740 // Make list of category columns to exclude from projection data
741 std::unique_ptr<RooArgSet> idxCompSliceSet( idxCatClone->getObservables(frame->getNormVars()) );
742
743 // If we don't project over the index, just do the regular plotOn
744 if (!projIndex) {
745
746 coutI(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") plot on " << frame->getPlotVar()->GetName()
747 << " represents a slice in the index category (" << _indexCat.arg().GetName() << ")" << endl ;
748
749 // Reduce projData: take out fitCat (component) columns and entries that don't match selected slice
750 // Construct cut string to only select projection data event that match the current slice
751
752 // Make cut string to exclude rows from projection data
753 TString cutString ;
754 bool first(true) ;
755 for (const auto arg : *idxCompSliceSet) {
756 auto idxComp = static_cast<RooCategory*>(arg);
757 RooAbsArg* slicedComponent = nullptr;
758 if (sliceSet && (slicedComponent = sliceSet->find(*idxComp)) != nullptr) {
759 auto theCat = static_cast<const RooAbsCategory*>(slicedComponent);
760 idxComp->setIndex(theCat->getCurrentIndex(), false);
761 }
762
763 if (!first) {
764 cutString.Append("&&") ;
765 } else {
766 first=false ;
767 }
768 cutString.Append(Form("%s==%d",idxComp->GetName(),idxComp->getCurrentIndex())) ;
769 }
770
771 // Make temporary projData without RooSim index category components
772 RooArgSet projDataVars(*projData->get()) ;
773 projDataVars.remove(*idxCompSliceSet,true,true) ;
774
775 std::unique_ptr<RooAbsData> projDataTmp( const_cast<RooAbsData*>(projData)->reduce(projDataVars,cutString) );
776
777 // Override normalization and projection dataset
778 RooCmdArg tmp1 = RooFit::Normalization(scaleFactor*wTable->getFrac(idxCatClone->getCurrentLabel()),stype) ;
779 RooCmdArg tmp2 = RooFit::ProjWData(*projDataSet,*projDataTmp) ;
780
781 // WVE -- do not adjust normalization for asymmetry plots
782 RooLinkedList cmdList2(cmdList) ;
783 if (!cmdList.find("Asymmetry")) {
784 cmdList2.Add(&tmp1) ;
785 }
786 cmdList2.Add(&tmp2) ;
787
788 // Plot single component
789 RooPlot* retFrame = getPdf(idxCatClone->getCurrentLabel())->plotOn(frame,cmdList2);
790 return retFrame ;
791 }
792
793 // If we project over the index, plot using a temporary RooAddPdf
794 // using the weights from the data as coefficients
795
796 // Build the list of indexCat components that are sliced
797 idxCompSliceSet->remove(projectedVars,true,true) ;
798
799 // Make a new expression that is the weighted sum of requested components
800 RooArgList pdfCompList ;
801 RooArgList wgtCompList ;
802//RooAbsPdf* pdf ;
803 double sumWeight(0) ;
804 for(auto * proxy : static_range_cast<RooRealProxy*>(_pdfProxyList)) {
805
806 idxCatClone->setLabel(proxy->name()) ;
807
808 // Determine if this component is the current slice (if we slice)
809 bool skip(false) ;
810 for (const auto idxSliceCompArg : *idxCompSliceSet) {
811 const auto idxSliceComp = static_cast<RooAbsCategory*>(idxSliceCompArg);
812 RooAbsCategory* idxComp = (RooAbsCategory*) idxCloneSet.find(idxSliceComp->GetName()) ;
813 if (idxComp->getCurrentIndex()!=idxSliceComp->getCurrentIndex()) {
814 skip=true ;
815 break ;
816 }
817 }
818 if (skip) continue ;
819
820 // Instantiate a RRV holding this pdfs weight fraction
821 wgtCompList.addOwned(std::make_unique<RooRealVar>(proxy->name(),"coef",wTable->getFrac(proxy->name())));
822 sumWeight += wTable->getFrac(proxy->name()) ;
823
824 // Add the PDF to list list
825 pdfCompList.add(proxy->arg()) ;
826 }
827
828 TString plotVarName(GetName()) ;
829 RooAddPdf plotVar{plotVarName,"weighted sum of RS components",pdfCompList,wgtCompList};
830
831 // Fix appropriate coefficient normalization in plot function
832 if (!_plotCoefNormSet.empty()) {
833 plotVar.fixAddCoefNormalization(_plotCoefNormSet) ;
834 }
835
836 std::unique_ptr<RooAbsData> projDataTmp;
837 RooArgSet projSetTmp ;
838 if (projData) {
839
840 // Construct cut string to only select projection data event that match the current slice
841 TString cutString ;
842 if (!idxCompSliceSet->empty()) {
843 bool first(true) ;
844 for (const auto idxSliceCompArg : *idxCompSliceSet) {
845 const auto idxSliceComp = static_cast<RooAbsCategory*>(idxSliceCompArg);
846 if (!first) {
847 cutString.Append("&&") ;
848 } else {
849 first=false ;
850 }
851 cutString.Append(Form("%s==%d",idxSliceComp->GetName(),idxSliceComp->getCurrentIndex())) ;
852 }
853 }
854
855 // Make temporary projData without RooSim index category components
856 RooArgSet projDataVars(*projData->get()) ;
857 RooArgSet idxCatServers;
858 _indexCat.arg().getObservables(frame->getNormVars(), idxCatServers) ;
859
860 projDataVars.remove(idxCatServers,true,true) ;
861
862 if (!idxCompSliceSet->empty()) {
863 projDataTmp.reset( const_cast<RooAbsData*>(projData)->reduce(projDataVars,cutString) );
864 } else {
865 projDataTmp.reset( const_cast<RooAbsData*>(projData)->reduce(projDataVars) );
866 }
867
868
869
870 if (projSet) {
871 projSetTmp.add(*projSet) ;
872 projSetTmp.remove(idxCatServers,true,true);
873 }
874 }
875
876
877 if (_indexCat.arg().isDerived() && !idxCompSliceSet->empty()) {
878 coutI(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") plot on " << frame->getPlotVar()->GetName()
879 << " represents a slice in index category components " << *idxCompSliceSet << endl ;
880
881 RooArgSet idxCompProjSet;
882 _indexCat.arg().getObservables(frame->getNormVars(), idxCompProjSet) ;
883 idxCompProjSet.remove(*idxCompSliceSet,true,true) ;
884 if (!idxCompProjSet.empty()) {
885 coutI(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") plot on " << frame->getPlotVar()->GetName()
886 << " averages with data index category components " << idxCompProjSet << endl ;
887 }
888 } else {
889 coutI(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") plot on " << frame->getPlotVar()->GetName()
890 << " averages with data index category (" << _indexCat.arg().GetName() << ")" << endl ;
891 }
892
893
894 // Override normalization and projection dataset
895 RooLinkedList cmdList2(cmdList) ;
896
897 RooCmdArg tmp1 = RooFit::Normalization(scaleFactor*sumWeight,stype) ;
898 RooCmdArg tmp2 = RooFit::ProjWData(*projDataSet,*projDataTmp) ;
899 // WVE -- do not adjust normalization for asymmetry plots
900 if (!cmdList.find("Asymmetry")) {
901 cmdList2.Add(&tmp1) ;
902 }
903 cmdList2.Add(&tmp2) ;
904
905 RooPlot* frame2 ;
906 if (!projSetTmp.empty()) {
907 // Plot temporary function
908 RooCmdArg tmp3 = RooFit::Project(projSetTmp) ;
909 cmdList2.Add(&tmp3) ;
910 frame2 = plotVar.plotOn(frame,cmdList2) ;
911 } else {
912 // Plot temporary function
913 frame2 = plotVar.plotOn(frame,cmdList2) ;
914 }
915
916 return frame2 ;
917}
918
919
920////////////////////////////////////////////////////////////////////////////////
921/// Interface function used by test statistics to freeze choice of observables
922/// for interpretation of fraction coefficients. Needed here because a RooSimultaneous
923/// works like a RooAddPdf when plotted
924
925void RooSimultaneous::selectNormalization(const RooArgSet* normSet, bool /*force*/)
926{
928 if (normSet) _plotCoefNormSet.add(*normSet) ;
929}
930
931
932////////////////////////////////////////////////////////////////////////////////
933/// Interface function used by test statistics to freeze choice of range
934/// for interpretation of fraction coefficients. Needed here because a RooSimultaneous
935/// works like a RooAddPdf when plotted
936
937void RooSimultaneous::selectNormalizationRange(const char* normRange2, bool /*force*/)
938{
939 _plotCoefNormRange = RooNameReg::ptr(normRange2) ;
940}
941
942
943
944
945////////////////////////////////////////////////////////////////////////////////
946
948 const RooArgSet* auxProto, bool verbose, bool autoBinned, const char* binnedTag) const
949{
950 const char* idxCatName = _indexCat.arg().GetName() ;
951
952 if (vars.find(idxCatName) && prototype==0
953 && (auxProto==0 || auxProto->empty())
954 && (autoBinned || (binnedTag && strlen(binnedTag)))) {
955
956 // Return special generator config that can also do binned generation for selected states
957 return new RooSimSplitGenContext(*this,vars,verbose,autoBinned,binnedTag) ;
958
959 } else {
960
961 // Return regular generator config ;
962 return genContext(vars,prototype,auxProto,verbose) ;
963 }
964}
965
966
967
968////////////////////////////////////////////////////////////////////////////////
969/// Return specialized generator context for simultaneous p.d.f.s
970
972 const RooArgSet* auxProto, bool verbose) const
973{
974 RooArgSet allVars{vars};
975 if(prototype) allVars.add(*prototype->get());
976
977 RooArgSet catsAmongAllVars;
978 allVars.selectCommon(flattenedCatList(), catsAmongAllVars);
979
980 // Not generating index cat: return context for pdf associated with present index state
981 if(catsAmongAllVars.empty()) {
982 auto* proxy = static_cast<RooRealProxy*>(_pdfProxyList.FindObject(_indexCat->getCurrentLabel()));
983 if (!proxy) {
984 coutE(InputArguments) << "RooSimultaneous::genContext(" << GetName()
985 << ") ERROR: no PDF associated with current state ("
986 << _indexCat.arg().GetName() << "=" << _indexCat.arg().getCurrentLabel() << ")" << endl ;
987 return nullptr;
988 }
989 return static_cast<RooAbsPdf*>(proxy->absArg())->genContext(vars,prototype,auxProto,verbose) ;
990 }
991
992 RooArgSet catsAmongProtoVars;
993 if(prototype) {
994 prototype->get()->selectCommon(flattenedCatList(), catsAmongProtoVars);
995
996 if(!catsAmongProtoVars.empty() && catsAmongProtoVars.size() != flattenedCatList().size()) {
997 // Abort if we have only part of the servers
998 coutE(Plotting) << "RooSimultaneous::genContext: ERROR: prototype must include either all "
999 << " components of the RooSimultaneous index category or none " << std::endl;
1000 return nullptr;
1001 }
1002 }
1003
1004 return new RooSimGenContext(*this,vars,prototype,auxProto,verbose) ;
1005}
1006
1007
1008
1009
1010////////////////////////////////////////////////////////////////////////////////
1011
1013 const RooArgSet* nset,
1014 double scaleFactor,
1015 bool correctForBinVolume,
1016 bool showProgress) const
1017{
1018 if (RooAbsReal::fillDataHist (hist, nset, scaleFactor,
1019 correctForBinVolume, showProgress) == 0)
1020 return 0;
1021
1022 const double sum = hist->sumEntries();
1023 if (sum != 0) {
1024 for (int i=0 ; i<hist->numEntries() ; i++) {
1025 hist->set(i, hist->weight(i) / sum, 0.);
1026 }
1027 }
1028
1029 return hist;
1030}
1031
1032
1033
1034
1035////////////////////////////////////////////////////////////////////////////////
1036/// Special generator interface for generation of 'global observables' -- for RooStats tools
1037
1039{
1040 // Make set with clone of variables (placeholder for output)
1041 RooArgSet globClone;
1042 whatVars.snapshot(globClone);
1043
1044 RooDataSet* data = new RooDataSet("gensimglobal","gensimglobal",whatVars) ;
1045
1046 for (Int_t i=0 ; i<nEvents ; i++) {
1047 for (const auto& nameIdx : indexCat()) {
1048
1049 // Get pdf associated with state from simpdf
1050 RooAbsPdf* pdftmp = getPdf(nameIdx.first);
1051
1052 // Generate only global variables defined by the pdf associated with this state
1053 RooArgSet globtmp;
1054 pdftmp->getObservables(&whatVars, globtmp) ;
1055 std::unique_ptr<RooDataSet> tmp{pdftmp->generate(globtmp,1)};
1056
1057 // Transfer values to output placeholder
1058 globClone.assign(*tmp->get(0)) ;
1059 }
1060 data->add(globClone) ;
1061 }
1062
1063 return data ;
1064}
1065
1066
1067/// Wraps the components of this RooSimultaneous in RooBinSamplingPdfs.
1068/// \param[in] data The dataset to be used in the eventual fit, used to figure
1069/// out the observables and whether the dataset is binned.
1070/// \param[in] precision Precision argument for all created RooBinSamplingPdfs.
1072
1073 if (precision < 0.) return;
1074
1075 RooArgSet newSamplingPdfs;
1076
1077 for (auto const &item : this->indexCat()) {
1078
1079 auto const &catName = item.first;
1080 auto &pdf = *this->getPdf(catName);
1081
1082 if (auto newSamplingPdf = RooBinSamplingPdf::create(pdf, data, precision)) {
1083 // Set the "ORIGNAME" attribute the indicate to
1084 // RooAbsArg::redirectServers() wich pdf should be replaced by this
1085 // RooBinSamplingPdf in the RooSimultaneous.
1086 newSamplingPdf->setAttribute(
1087 (std::string("ORIGNAME:") + pdf.GetName()).c_str());
1088 newSamplingPdfs.addOwned(std::move(newSamplingPdf));
1089 }
1090 }
1091
1092 this->redirectServers(newSamplingPdfs, false, true);
1093 this->addOwnedComponents(std::move(newSamplingPdfs));
1094}
1095
1096
1097/// Wraps the components of this RooSimultaneous in RooBinSamplingPdfs, with a
1098/// different precision parameter for each component.
1099/// \param[in] data The dataset to be used in the eventual fit, used to figure
1100/// out the observables and whether the dataset is binned.
1101/// \param[in] precisions The map that gives the precision argument for each
1102/// component in the RooSimultaneous. The keys are the pdf names. If
1103/// there is no value for a given component, it will not use the bin
1104/// integration. Otherwise, the value has the same meaning than in
1105/// the IntegrateBins() command argument for RooAbsPdf::fitTo().
1106/// \param[in] useCategoryNames If this flag is set, the category names will be
1107/// used to look up the precision in the precisions map instead of
1108/// the pdf names.
1110 std::map<std::string, double> const& precisions,
1111 bool useCategoryNames /*=false*/) {
1112
1113 constexpr double defaultPrecision = -1.;
1114
1115 RooArgSet newSamplingPdfs;
1116
1117 for (auto const &item : this->indexCat()) {
1118
1119 auto const &catName = item.first;
1120 auto &pdf = *this->getPdf(catName);
1121 std::string pdfName = pdf.GetName();
1122
1123 auto found = precisions.find(useCategoryNames ? catName : pdfName);
1124 const double precision =
1125 found != precisions.end() ? found->second : defaultPrecision;
1126 if (precision < 0.)
1127 continue;
1128
1129 if (auto newSamplingPdf = RooBinSamplingPdf::create(pdf, data, precision)) {
1130 // Set the "ORIGNAME" attribute the indicate to
1131 // RooAbsArg::redirectServers() wich pdf should be replaced by this
1132 // RooBinSamplingPdf in the RooSimultaneous.
1133 newSamplingPdf->setAttribute(
1134 (std::string("ORIGNAME:") + pdf.GetName()).c_str());
1135 newSamplingPdfs.addOwned(std::move(newSamplingPdf));
1136 }
1137 }
1138
1139 this->redirectServers(newSamplingPdfs, false, true);
1140 this->addOwnedComponents(std::move(newSamplingPdfs));
1141}
1142
1143/// Internal utility function to get a list of all category components for this
1144/// RooSimultaneous. The output contains only the index category if it is a
1145/// RooCategory, or the list of all category components if it is a
1146/// RooSuperCategory.
1148{
1149 // Note that the index category of a RooSimultaneous can only be of type
1150 // RooCategory or RooSuperCategory, because these are the only classes that
1151 // inherit from RooAbsCategoryLValue.
1152 if (auto superCat = dynamic_cast<RooSuperCategory const*>(&_indexCat.arg())) {
1153 return superCat->inputCatList();
1154 }
1155
1156 if(!_indexCatSet) {
1157 _indexCatSet = std::make_unique<RooArgSet>(_indexCat.arg());
1158 }
1159 return *_indexCatSet;
1160}
1161
1162namespace {
1163
1164void prefixArgs(RooAbsArg *arg, std::string const &prefix, RooArgSet const &normSet)
1165{
1166 if (!arg->getStringAttribute("__prefix__")) {
1167 arg->SetName((prefix + arg->GetName()).c_str());
1168 arg->setStringAttribute("__prefix__", prefix.c_str());
1169 }
1170 for (RooAbsArg *server : arg->servers()) {
1171 if (server->isFundamental() && normSet.find(*server)) {
1172 prefixArgs(server, prefix, normSet);
1173 server->setAttribute("__obs__");
1174 } else if (!server->isFundamental()) {
1175 prefixArgs(server, prefix, normSet);
1176 }
1177 }
1178}
1179
1180} // namespace
1181
1182std::unique_ptr<RooAbsArg>
1184{
1185 std::unique_ptr<RooSimultaneous> newSimPdf{static_cast<RooSimultaneous *>(this->Clone())};
1186
1187 const char *rangeName = this->getStringAttribute("RangeName");
1188 bool splitRange = this->getAttribute("SplitRange");
1189
1190 RooArgSet newPdfs;
1191 std::vector<std::string> catNames;
1192
1193 for (auto *proxy : static_range_cast<RooRealProxy *>(newSimPdf->_pdfProxyList)) {
1194 catNames.emplace_back(proxy->GetName());
1195 std::string const &catName = catNames.back();
1196 const std::string prefix = "_" + catName + "_";
1197
1198 const std::string origname = proxy->arg().GetName();
1199
1200 auto pdfClone = RooHelpers::cloneTreeWithSameParameters(static_cast<RooAbsPdf const &>(proxy->arg()), &normSet);
1201
1202 prefixArgs(pdfClone.get(), prefix, normSet);
1203
1204 auto binnedInfo = RooHelpers::getBinnedL(*pdfClone);
1205
1206 RooAbsPdf &pdf = binnedInfo.binnedPdf ? *binnedInfo.binnedPdf : *pdfClone;
1207
1208 if (binnedInfo.isBinnedL) {
1209 pdf.setAttribute("BinnedLikelihoodActive");
1210 }
1211
1212 std::unique_ptr<RooArgSet> pdfNormSet(
1213 static_cast<RooArgSet *>(std::unique_ptr<RooArgSet>(pdf.getVariables())->selectByAttrib("__obs__", true)));
1214
1215 if (rangeName) {
1216 pdf.setNormRange(RooHelpers::getRangeNameForSimComponent(rangeName, splitRange, catName).c_str());
1217 }
1218
1219 auto *pdfFinal = RooFit::Detail::CompileContext{*pdfNormSet}.compile(pdf, *newSimPdf, *pdfNormSet);
1220 pdfFinal->fixAddCoefNormalization(*pdfNormSet, false);
1221
1222 pdfClone->SetName((std::string("_") + pdfClone->GetName()).c_str());
1223 pdfFinal->addOwnedComponents(std::move(pdfClone));
1224
1225 pdfFinal->setAttribute(("ORIGNAME:" + origname).c_str());
1226 newPdfs.add(*pdfFinal);
1227
1228 // We will remove the old pdf server because we will fill the new ones by
1229 // hand via the creation of new proxies.
1230 newSimPdf->removeServer(const_cast<RooAbsReal &>(proxy->arg()), true);
1231 }
1232
1233 // Replace pdfs with compiled pdfs. Don't use RooAbsArg::redirectServers()
1234 // here, because it doesn't support replacing two servers with the same name
1235 // (it can happen in a RooSimultaneous that two pdfs have the same name).
1236
1237 // First delete old proxies (we have already removed the servers before).
1238 newSimPdf->_pdfProxyList.Delete();
1239
1240 // Recreate the _pdfProxyList with the compiled pdfs
1241 for (std::size_t i = 0; i < newPdfs.size(); ++i) {
1242 const char *label = catNames[i].c_str();
1243 newSimPdf->_pdfProxyList.Add(
1244 new RooRealProxy(label, label, newSimPdf.get(), *static_cast<RooAbsReal *>(newPdfs[i])));
1245 }
1246
1247 ctx.compileServers(*newSimPdf, normSet); // to trigger compiling also the index category
1248
1249 return newSimPdf;
1250}
#define coutI(a)
#define oocoutW(o, a)
#define oocxcoutD(o, a)
#define oocoutI(o, a)
#define coutE(a)
RooTemplateProxy< RooAbsReal > RooRealProxy
Compatibility typedef replacing the old RooRealProxy class.
#define ClassImp(name)
Definition Rtypes.h:377
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:2468
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition RooAbsArg.h:78
bool dependsOn(const RooAbsCollection &serverList, const RooAbsArg *ignoreArg=nullptr, bool valueOnly=false) const
Test whether we depend on (ie, are served by) any object in the specified collection.
friend class RooDataSet
Definition RooAbsArg.h:669
void setStringAttribute(const Text_t *key, const Text_t *value)
Associate string 'value' to this object under key 'key'.
RooFit::OwningPtr< RooArgSet > getObservables(const RooArgSet &set, bool valueOnly=true) const
Given a set of possible observables, return the observables that this PDF depends on.
void SetName(const char *name) override
Set the name of the TNamed.
bool addOwnedComponents(const RooAbsCollection &comps)
Take ownership of the contents of 'comps'.
const Text_t * getStringAttribute(const Text_t *key) const
Get string attribute mapped under key 'key'.
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.
virtual bool isDerived() const
Does value or shape of this arg depend on any other arg?
Definition RooAbsArg.h:98
bool getAttribute(const Text_t *name) const
Check if a named attribute is set. By default, all attributes are unset.
RooFit::OwningPtr< RooArgSet > getVariables(bool stripDisconnected=true) const
Return RooArgSet with all variables (tree leaf nodes of expresssion tree)
void setAttribute(const Text_t *name, bool value=true)
Set (default) or clear a named boolean attribute of this object.
TObject * Clone(const char *newname=nullptr) const override
Make a clone of an object using the Streamer facility.
Definition RooAbsArg.h:90
virtual bool isFundamental() const
Is this object a fundamental type that can be added to a dataset? Fundamental-type subclasses overrid...
Definition RooAbsArg.h:251
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 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.
const char * GetName() const override
Returns name of object.
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
TObject * Clone(const char *newname=nullptr) const override
Make a clone of an object using the Streamer facility.
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
virtual bool addOwned(RooAbsArg &var, bool silent=false)
Add an argument and transfer the ownership to the collection.
bool selectCommon(const RooAbsCollection &refColl, RooAbsCollection &outColl) const
Create a subset of the current collection, consisting only of those elements that are contained as we...
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:58
virtual const RooArgSet * get() const
Definition RooAbsData.h:102
virtual Roo1DTable * table(const RooArgSet &catSet, const char *cuts="", const char *opts="") const
Construct table for product of categories in catSet.
virtual Int_t numEntries() const
Return number of entries in dataset, i.e., count unweighted entries.
RooAbsGenContext is the abstract base class for generator contexts of RooAbsPdf objects.
void setNormRange(const char *rangeName)
bool mustBeExtended() const
If true PDF must provide extended likelihood term.
Definition RooAbsPdf.h:275
RooArgSet const * _normSet
Normalization integral (owned by _normMgr)
Definition RooAbsPdf.h:368
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:271
@ CanBeExtended
Definition RooAbsPdf.h:265
@ MustBeExtended
Definition RooAbsPdf.h:265
@ CanNotBeExtended
Definition RooAbsPdf.h:265
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:61
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:105
bool plotSanityChecks(RooPlot *frame) const
Utility function for plotOn(), perform general sanity check on frame to ensure safe plotting operatio...
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:33
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:55
RooArgSet * snapshot(bool deepCopy=true) const
Use RooAbsCollection::snapshot(), but return as RooArgSet.
Definition RooArgSet.h:178
static std::unique_ptr< RooAbsPdf > create(RooAbsPdf &pdf, RooAbsData const &data, double precision)
Creates a wrapping RooBinSamplingPdf if appropriate.
Int_t setObj(const RooArgSet *nset, T *obj, const TNamed *isetRangeName=nullptr)
Setter 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.
T * getObj(const RooArgSet *nset, Int_t *sterileIndex=nullptr, const TNamed *isetRangeName=nullptr)
Getter 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.
bool defineObject(const char *name, const char *argName, Int_t setNum, const TObject *obj=nullptr, bool isArray=false)
Define TObject property name 'name' mapped to object in slot 'setNum' in RooCmdArg with name argName ...
double getDouble(const char *name, double defaultValue=0.0)
Return double property registered with name 'name'.
RooArgSet * getSet(const char *name, RooArgSet *set=nullptr)
Return RooArgSet property registered with name 'name'.
void defineMutex(const char *head, Args_t &&... tail)
Define arguments where any pair is mutually exclusive.
bool process(const RooCmdArg &arg)
Process given RooCmdArg.
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 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...
bool defineString(const char *name, const char *argName, Int_t stringNum, const char *defValue="", bool appendMode=false)
Define double property name 'name' mapped to double in slot 'stringNum' in RooCmdArg with name argNam...
bool ok(bool verbose) const
Return true of parsing was successful.
TObject * getObject(const char *name, TObject *obj=nullptr)
Return TObject property registered with name 'name'.
const char * getString(const char *name, const char *defaultValue="", bool convEmptyToNull=false)
Return string property registered with name 'name'.
void allowUndefined(bool flag=true)
If flag is true the processing of unrecognized RooCmdArgs is not considered an error.
bool defineSet(const char *name, const char *argName, Int_t setNum, const RooArgSet *set=nullptr)
Define TObject property name 'name' mapped to object in slot 'setNum' in RooCmdArg with name argName ...
bool defineDouble(const char *name, const char *argName, Int_t doubleNum, double defValue=0.0)
Define double property name 'name' mapped to double in slot 'doubleNum' in RooCmdArg with name argNam...
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:39
double weight(std::size_t i) const
Return weight of i-th bin.
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 sumEntries() const override
Sum the weights of all bins.
RooDataSet is a container class to hold unbinned data.
Definition RooDataSet.h:57
const RooArgSet * get(Int_t index) const override
Return RooArgSet with coordinates of event 'index'.
void compileServers(RooAbsArg &arg, RooArgSet const &normSet)
T * compile(T &arg, RooAbsArg &owner, RooArgSet const &normSet)
RooLinkedList is an collection class for internal use, storing a collection of RooAbsArg pointers in ...
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:43
void remove(const char *name=nullptr, bool deleteToo=true)
Remove object with given name, or last object added if no name is given.
Definition RooPlot.cxx:883
const RooArgSet * getNormVars() const
Definition RooPlot.h:146
RooAbsRealLValue * getPlotVar() const
Definition RooPlot.h:137
RooSimultaneous facilitates simultaneous fitting of multiple PDFs to subsets of a given dataset.
void selectNormalization(const RooArgSet *depSet=nullptr, bool force=false) override
Interface function used by test statistics to freeze choice of observables for interpretation of frac...
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.
void selectNormalizationRange(const char *rangeName=nullptr, bool force=false) override
Interface function used by test statistics to freeze choice of range for interpretation of fraction c...
TList _pdfProxyList
List of PDF proxies (named after applicable category state)
RooObjCacheManager _partIntMgr
! Component normalization manager
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.
~RooSimultaneous() override
Destructor.
RooArgSet const & flattenedCatList() const
Internal utility function to get a list of all category components for this RooSimultaneous.
ExtendMode extendMode() const override
Examine the pdf components and check if one of them can be extended or must be extended It is enough ...
RooCategoryProxy _indexCat
Index category.
double analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=nullptr) const override
Return analytical integration defined by given code.
Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &numVars, const RooArgSet *normSet, const char *rangeName=nullptr) const override
Forward determination of analytical integration capabilities to component p.d.f.s A unique code is as...
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
RooAbsGenContext * autoGenContext(const RooArgSet &vars, const RooDataSet *prototype=nullptr, const RooArgSet *auxProto=nullptr, bool verbose=false, bool autoBinned=true, const char *binnedTag="") const override
RooAbsPdf * getPdf(RooStringView catName) const
Return the p.d.f associated with the given index category name.
std::unique_ptr< RooAbsArg > compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext &ctx) const override
static std::unique_ptr< RooSimultaneous::InitializationOutput > initialize(std::string const &name, RooAbsCategoryLValue &inIndexCat, std::map< std::string, RooAbsPdf * > const &pdfMap)
void wrapPdfsInBinSamplingPdfs(RooAbsData const &data, double precision)
Wraps the components of this RooSimultaneous in RooBinSamplingPdfs.
const TNamed * _plotCoefNormRange
RooAbsGenContext * genContext(const RooArgSet &vars, const RooDataSet *prototype=nullptr, const RooArgSet *auxProto=nullptr, bool verbose=false) const override
Return specialized generator context for simultaneous p.d.f.s.
std::unique_ptr< RooArgSet > _indexCatSet
! Index category wrapped in a RooArgSet if needed internally
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
The RooStringView is a wrapper around a C-syle string that can also be constructed from a std::string...
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.
const T & arg() const
Return reference to object held in proxy.
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
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:41
virtual Int_t IndexOf(const TObject *obj) const
Return index of object in collection.
Basic string class.
Definition TString.h:139
TString & Append(const char *cs)
Definition TString.h:576
RooCmdArg ProjWData(const RooAbsData &projData, bool binData=false)
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.
std::unique_ptr< T > cloneTreeWithSameParameters(T const &arg, RooArgSet const *observables=nullptr)
Clone RooAbsArg object and reattach to original parameters.
Definition RooHelpers.h:149
BinnedLOutput getBinnedL(RooAbsPdf &pdf)
std::string getRangeNameForSimComponent(std::string const &rangeName, bool splitRange, std::string const &catName)
Definition first.py:1
std::unique_ptr< RooArgSet > subIndexComps
const RooAbsCategoryLValue * subIndex
Internal struct used for initialization.
std::vector< RooAbsPdf const * > finalPdfs
std::vector< std::string > finalCatLabels
void addPdf(const RooAbsPdf &pdf, std::string const &catLabel)
static uint64_t sum(uint64_t i)
Definition Factory.cxx:2345