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
22Facilitates simultaneous fitting of multiple PDFs to subsets of a given dataset.
23The class takes an index category, which is used as a selector
24for PDFs, and a list of PDFs, each associated
25with a state of the index category. RooSimultaneous always returns
26the value of the PDF that is associated with the current value
27of the index category.
28
29Extended likelihood fitting is supported if all components support
30extended likelihood mode. The expected number of events by a RooSimultaneous
31is that of the component p.d.f. selected by the index category.
32
33The index category can be accessed using indexCategory().
34
35###Generating events
36When generating events from a RooSimultaneous, the index category has to be added to
37the dataset. Further, the PDF needs to know the relative probabilities of each category, i.e.,
38how many events are in which category. This can be achieved in two ways:
39- Generating with proto data that have category entries: An event from the same category as
40in the proto data is created for each event in the proto data.
41See RooAbsPdf::generate(const RooArgSet&,const RooDataSet&,Int_t,bool,bool,bool) const.
42- No proto data: A category is chosen randomly.
43\note This requires that the PDFs building the simultaneous are extended. In this way,
44the relative probability of each category can be calculated from the number of events
45in each category.
46**/
47
48#include "RooSimultaneous.h"
49
50#include "Roo1DTable.h"
52#include "RooAbsData.h"
53#include "RooAddPdf.h"
54#include "RooArgSet.h"
55#include "RooBinSamplingPdf.h"
56#include "RooCategory.h"
57#include "RooCmdConfig.h"
58#include "RooDataHist.h"
59#include "RooDataSet.h"
60#include "RooGlobalFunc.h"
61#include "RooMsgService.h"
62#include "RooNameReg.h"
63#include "RooPlot.h"
64#include "RooRandom.h"
65#include "RooRealVar.h"
66#include "RooSimGenContext.h"
68#include "RooSuperCategory.h"
69
70#include "RooFitImplHelpers.h"
71
72#include <ROOT/StringUtils.hxx>
73
74#include <iostream>
75
76namespace {
77
78std::map<std::string, RooAbsPdf *> createPdfMap(const RooArgList &inPdfList, RooAbsCategoryLValue &inIndexCat)
79{
80 std::map<std::string, RooAbsPdf *> pdfMap;
82 for (unsigned int i = 0; i < inPdfList.size(); ++i) {
83 auto pdf = static_cast<RooAbsPdf *>(&inPdfList[i]);
84 const auto &nameIdx = (*indexCatIt++);
85 pdfMap[nameIdx.first] = pdf;
86 }
87 return pdfMap;
88}
89
90} // namespace
91
93
95{
96 finalPdfs.push_back(&pdf);
97 finalCatLabels.emplace_back(catLabel);
98}
99
100using std::string;
101
102
103
104
105////////////////////////////////////////////////////////////////////////////////
106/// Constructor with index category. PDFs associated with indexCat
107/// states can be added after construction with the addPdf() function.
108///
109/// RooSimultaneous can function without having a PDF associated
110/// with every single state. The normalization in such cases is taken
111/// from the number of registered PDFs, but getVal() will assert if
112/// when called for an unregistered index state.
113
114RooSimultaneous::RooSimultaneous(const char *name, const char *title,
116 RooSimultaneous{name, title, std::map<std::string, RooAbsPdf*>{}, inIndexCat}
117{
118}
119
120
121////////////////////////////////////////////////////////////////////////////////
122/// Constructor from index category and full list of PDFs.
123/// In this constructor form, a PDF must be supplied for each indexCat state
124/// to avoid ambiguities. The PDFs are associated with the states of the
125/// index category as they appear when iterating through the category states
126/// with RooAbsCategory::begin() and RooAbsCategory::end(). This usually means
127/// they are associated by ascending index numbers.
128///
129/// PDFs may not overlap (i.e. share any variables) with the index category (function)
130
131RooSimultaneous::RooSimultaneous(const char *name, const char *title,
134{
135 if (inPdfList.size() != inIndexCat.size()) {
136 std::stringstream errMsg;
137 errMsg << "RooSimultaneous::ctor(" << GetName()
138 << " ERROR: Number PDF list entries must match number of index category states, no PDFs added";
139 coutE(InputArguments) << errMsg.str() << std::endl;
140 throw std::invalid_argument(errMsg.str());
141 }
142}
143
144
145////////////////////////////////////////////////////////////////////////////////
146
147RooSimultaneous::RooSimultaneous(const char *name, const char *title, std::map<string, RooAbsPdf *> pdfMap,
149 : RooSimultaneous(name, title, std::move(*initialize(name ? name : "", inIndexCat, pdfMap)))
150{
151}
152
153/// For internal use in RooFit.
154RooSimultaneous::RooSimultaneous(const char *name, const char *title,
157 : RooSimultaneous(name, title, RooFit::Detail::flatMapToStdMap(pdfMap), inIndexCat)
158{
159}
160
162 : RooAbsPdf(name, title),
163 _plotCoefNormSet("!plotCoefNormSet", "plotCoefNormSet", this, false, false),
164 _partIntMgr(this, 10),
165 _indexCat("indexCat", "Index category", this, *initInfo.indexCat)
166{
167 for (std::size_t i = 0; i < initInfo.finalPdfs.size(); ++i) {
168 addPdf(*initInfo.finalPdfs[i], initInfo.finalCatLabels[i].c_str());
169 }
170
171 // Take ownership of eventual super category
172 if (initInfo.superIndex) {
173 addOwnedComponents(std::move(initInfo.superIndex));
174 }
175}
176
177/// \cond ROOFIT_INTERNAL
178
179// This class cannot be locally defined in initialize as it cannot be
180// used as a template argument in that case
181namespace RooSimultaneousAux {
182 struct CompInfo {
183 RooAbsPdf* pdf ;
186 std::unique_ptr<RooArgSet> subIndexComps;
187 } ;
188}
189
190/// \endcond
191
192std::unique_ptr<RooSimultaneous::InitializationOutput>
194 std::map<std::string, RooAbsPdf *> const& pdfMap)
195
196{
197 auto out = std::make_unique<RooSimultaneous::InitializationOutput>();
198 out->indexCat = &inIndexCat;
199
200 // First see if there are any RooSimultaneous input components
201 bool simComps(false) ;
202 for (auto const& item : pdfMap) {
203 if (dynamic_cast<RooSimultaneous*>(item.second)) {
204 simComps = true ;
205 break ;
206 }
207 }
208
209 // If there are no simultaneous component p.d.f. do simple processing through addPdf()
210 if (!simComps) {
211 for (auto const& item : pdfMap) {
212 out->addPdf(*item.second,item.first);
213 }
214 return out;
215 }
216
217 std::string msgPrefix = "RooSimultaneous::initialize(" + name + ") ";
218
219 // Issue info message that we are about to do some rearranging
220 oocoutI(nullptr, InputArguments) << msgPrefix << "INFO: one or more input component of simultaneous p.d.f.s are"
221 << " simultaneous p.d.f.s themselves, rewriting composite expressions as one-level simultaneous p.d.f. in terms of"
222 << " final constituents and extended index category" << std::endl;
223
224
226 std::map<string,RooSimultaneousAux::CompInfo> compMap ;
227 for (auto const& item : pdfMap) {
228 RooSimultaneousAux::CompInfo ci ;
229 ci.pdf = item.second ;
230 RooSimultaneous* simComp = dynamic_cast<RooSimultaneous*>(item.second) ;
231 if (simComp) {
232 ci.simPdf = simComp ;
233 ci.subIndex = &simComp->indexCat() ;
234 ci.subIndexComps = simComp->indexCat().isFundamental()
235 ? std::make_unique<RooArgSet>(simComp->indexCat())
236 : std::unique_ptr<RooArgSet>(simComp->indexCat().getVariables());
237 allAuxCats.add(*ci.subIndexComps,true) ;
238 } else {
239 ci.simPdf = nullptr;
240 ci.subIndex = nullptr;
241 }
242 compMap[item.first] = std::move(ci);
243 }
244
245 // Construct the 'superIndex' from the nominal index category and all auxiliary components
246 RooArgSet allCats(inIndexCat) ;
247 allCats.add(allAuxCats) ;
248 std::string siname = name + "_index";
249 out->superIndex = std::make_unique<RooSuperCategory>(siname.c_str(),siname.c_str(),allCats) ;
250 auto *superIndex = out->superIndex.get();
251 out->indexCat = superIndex;
252
253 // Now process each of original pdf/state map entries
254 for (auto const& citem : compMap) {
255
257 if (citem.second.subIndexComps) {
258 repliCats.remove(*citem.second.subIndexComps) ;
259 }
260 inIndexCat.setLabel(citem.first.c_str()) ;
261
262 if (!citem.second.simPdf) {
263
264 // Entry is a plain p.d.f. assign it to every state permutation of the repliCats set
266
267 // Iterator over all states of repliSuperCat
268 for (const auto& nameIdx : repliSuperCat) {
269 // Set value
270 repliSuperCat.setLabel(nameIdx.first) ;
271 // Retrieve corresponding label of superIndex
272 string superLabel = superIndex->getCurrentLabel() ;
273 out->addPdf(*citem.second.pdf,superLabel);
274 oocxcoutD(static_cast<RooAbsArg*>(nullptr), InputArguments) << msgPrefix
275 << "assigning pdf " << citem.second.pdf->GetName() << " to super label " << superLabel << std::endl ;
276 }
277 } else {
278
279 // Entry is a simultaneous p.d.f
280
281 if (repliCats.empty()) {
282
283 // Case 1 -- No replication of components of RooSim component are required
284
285 for (const auto& type : *citem.second.subIndex) {
286 const_cast<RooAbsCategoryLValue*>(citem.second.subIndex)->setLabel(type.first.c_str());
287 string superLabel = superIndex->getCurrentLabel() ;
288 RooAbsPdf* compPdf = citem.second.simPdf->getPdf(type.first);
289 if (compPdf) {
290 out->addPdf(*compPdf,superLabel);
291 oocxcoutD(static_cast<RooAbsArg*>(nullptr), InputArguments) << msgPrefix
292 << "assigning pdf " << compPdf->GetName() << "(member of " << citem.second.pdf->GetName()
293 << ") to super label " << superLabel << std::endl ;
294 } else {
295 oocoutW(nullptr, InputArguments) << msgPrefix << "WARNING: No p.d.f. associated with label "
296 << type.second << " for component RooSimultaneous p.d.f " << citem.second.pdf->GetName()
297 << "which is associated with master index label " << citem.first << std::endl ;
298 }
299 }
300
301 } else {
302
303 // Case 2 -- Replication of components of RooSim component are required
304
305 // Make replication supercat
307
308 for (const auto& stype : *citem.second.subIndex) {
309 const_cast<RooAbsCategoryLValue*>(citem.second.subIndex)->setLabel(stype.first.c_str());
310
311 for (const auto& nameIdx : repliSuperCat) {
312 repliSuperCat.setLabel(nameIdx.first) ;
313 const string superLabel = superIndex->getCurrentLabel() ;
314 RooAbsPdf* compPdf = citem.second.simPdf->getPdf(stype.first);
315 if (compPdf) {
316 out->addPdf(*compPdf,superLabel);
317 oocxcoutD(static_cast<RooAbsArg*>(nullptr), InputArguments) << msgPrefix
318 << "assigning pdf " << compPdf->GetName() << "(member of " << citem.second.pdf->GetName()
319 << ") to super label " << superLabel << std::endl ;
320 } else {
321 oocoutW(nullptr, InputArguments) << msgPrefix << "WARNING: No p.d.f. associated with label "
322 << stype.second << " for component RooSimultaneous p.d.f " << citem.second.pdf->GetName()
323 << "which is associated with master index label " << citem.first << std::endl ;
324 }
325 }
326 }
327 }
328 }
329 }
330
331 return out;
332}
333
334
335////////////////////////////////////////////////////////////////////////////////
336/// Copy constructor
337
340 _plotCoefNormSet("!plotCoefNormSet",this,other._plotCoefNormSet),
341 _plotCoefNormRange(other._plotCoefNormRange),
342 _partIntMgr(other._partIntMgr,this),
343 _indexCat("indexCat",this,other._indexCat),
344 _numPdf(other._numPdf)
345{
346 // Copy proxy list
347 for(auto* proxy : static_range_cast<RooRealProxy*>(other._pdfProxyList)) {
348 _pdfProxyList.Add(new RooRealProxy(proxy->GetName(),this,*proxy)) ;
349 }
350}
351
352
353
354////////////////////////////////////////////////////////////////////////////////
355/// Destructor
356
361
362
363
364////////////////////////////////////////////////////////////////////////////////
365/// Return the p.d.f associated with the given index category name
366
368{
370 return proxy ? static_cast<RooAbsPdf*>(proxy->absArg()) : nullptr;
371}
372
373
374
375////////////////////////////////////////////////////////////////////////////////
376/// Associate given PDF with index category state label 'catLabel'.
377/// The name state must be already defined in the index category.
378///
379/// RooSimultaneous can function without having a PDF associated
380/// with every single state. The normalization in such cases is taken
381/// from the number of registered PDFs, but getVal() will fail if
382/// called for an unregistered index state.
383///
384/// PDFs may not overlap (i.e. share any variables) with the index category (function).
385/// \param[in] pdf PDF to be added.
386/// \param[in] catLabel Name of the category state to be associated to the PDF.
387/// \return `true` in case of failure.
388
389bool RooSimultaneous::addPdf(const RooAbsPdf& pdf, const char* catLabel)
390{
391 // PDFs cannot overlap with the index category
392 if (pdf.dependsOn(_indexCat.arg())) {
393 coutE(InputArguments) << "RooSimultaneous::addPdf(" << GetName() << "): PDF '" << pdf.GetName()
394 << "' overlaps with index category '" << _indexCat.arg().GetName() << "'."<< std::endl ;
395 return true ;
396 }
397
398 // Each index state can only have one PDF associated with it
400 coutE(InputArguments) << "RooSimultaneous::addPdf(" << GetName() << "): index state '"
401 << catLabel << "' has already an associated PDF." << std::endl ;
402 return true ;
403 }
404
405 const RooSimultaneous* simPdf = dynamic_cast<const RooSimultaneous*>(&pdf) ;
406 if (simPdf) {
407
408 coutE(InputArguments) << "RooSimultaneous::addPdf(" << GetName()
409 << ") ERROR: you cannot add a RooSimultaneous component to a RooSimultaneous using addPdf()."
410 << " Use the constructor with RooArgList if input p.d.f.s or the map<string,RooAbsPdf&> instead." << std::endl ;
411 return true ;
412
413 } else {
414
415 // Create a proxy named after the associated index state
416 TObject* proxy = new RooRealProxy(catLabel,catLabel,this,const_cast<RooAbsPdf&>(pdf));
418 _numPdf += 1 ;
419 }
420
421 return false ;
422}
423
424////////////////////////////////////////////////////////////////////////////////
425/// Examine the pdf components and check if one of them can be extended or must be extended.
426/// It is enough to have one component that can be extended or must be extended to return the flag in
427/// the total simultaneous pdf.
428
430{
431 bool anyCanExtend = false;
432
434 auto &pdf = static_cast<RooAbsPdf const&>(proxy->arg());
435 if (pdf.mustBeExtended())
436 return MustBeExtended;
437 anyCanExtend |= pdf.canBeExtended();
438 }
440}
441
442////////////////////////////////////////////////////////////////////////////////
443/// Return the current value:
444/// the value of the PDF associated with the current index category state
445
447{
448 // Retrieve the proxy by index name
450 if(!proxy) {
451 return 0;
452 }
453
454 double nEvtTot = 1.0;
455 double nEvtCat = 1.0;
456
457 // Calculate relative weighting factor for sim-pdfs of all extendable components
458 if (canBeExtended()) {
459
460 nEvtTot = 0;
461 nEvtCat = 0;
462
464 auto &pdf2 = static_cast<RooAbsPdf const &>(proxy2->arg());
465 if(!pdf2.canBeExtended()) {
466 // If one of the pdfs can't be expected, reset the normalization
467 // factor to one and break out of the loop.
468 nEvtTot = 1.0;
469 nEvtCat = 1.0;
470 break;
471 }
472 const double nEvt = pdf2.expectedEvents(_normSet);
473 nEvtTot += nEvt;
474 if (proxy == proxy2) {
475 // Matching by proxy by pointer rather than pdfs, because it's
476 // possible to have the same pdf used in different states.
477 nEvtCat += nEvt;
478 }
479 }
480 }
481 double catFrac = nEvtCat / nEvtTot;
482
483 // Return the selected PDF value, normalized by the relative number of
484 // expected events if applicable.
485 return *proxy * catFrac;
486}
487
488////////////////////////////////////////////////////////////////////////////////
489/// Return the number of expected events: If the index is in nset,
490/// then return the sum of the expected events of all components,
491/// otherwise return the number of expected events of the PDF
492/// associated with the current index category state
493
495{
496 if (nset->contains(_indexCat.arg())) {
497
498 double sum(0) ;
499
501 sum += (static_cast<RooAbsPdf*>(proxy->absArg()))->expectedEvents(nset) ;
502 }
503
504 return sum ;
505
506 } else {
507
508 // Retrieve the proxy by index name
510
511 //assert(proxy!=0) ;
512 if (proxy==nullptr) return 0 ;
513
514 // Return the selected PDF value, normalized by the number of index states
515 return (static_cast<RooAbsPdf*>(proxy->absArg()))->expectedEvents(nset);
516 }
517}
518
519
520
521////////////////////////////////////////////////////////////////////////////////
522/// Forward determination of analytical integration capabilities to component p.d.f.s
523/// A unique code is assigned to the combined integration capabilities of all associated
524/// p.d.f.s
525
527 const RooArgSet* normSet, const char* rangeName) const
528{
529 // Declare that we can analytically integrate all requested observables
530 analVars.add(allVars) ;
531
532 // Retrieve (or create) the required partial integral list
533 Int_t code ;
534
535 // Check if this configuration was created before
537 if (cache) {
538 code = _partIntMgr.lastIndex() ;
539 return code+1 ;
540 }
541 cache = new CacheElem ;
542
543 // Create the partial integral set for this request
545 cache->_partIntList.addOwned(std::unique_ptr<RooAbsReal>{proxy->arg().createIntegral(analVars,normSet,nullptr,rangeName)});
546 }
547
548 // Store the partial integral list and return the assigned code ;
550
551 return code+1 ;
552}
553
554
555
556////////////////////////////////////////////////////////////////////////////////
557/// Return analytical integration defined by given code
558
559double RooSimultaneous::analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* /*rangeName*/) const
560{
561 // No integration scenario
562 if (code==0) {
563 return getVal(normSet) ;
564 }
565
566 // Partial integration scenarios, rangeName already encoded in 'code'
567 CacheElem* cache = static_cast<CacheElem*>(_partIntMgr.getObjByIndex(code-1)) ;
568
571 return (static_cast<RooAbsReal*>(cache->_partIntList.at(idx)))->getVal(normSet) ;
572}
573
574
575
576
577
578
579////////////////////////////////////////////////////////////////////////////////
580/// Back-end for plotOn() implementation on RooSimultaneous which
581/// needs special handling because a RooSimultaneous PDF cannot
582/// project out its index category via integration. plotOn() will
583/// abort if this is requested without providing a projection dataset.
584
586{
587 // Sanity checks
588 if (plotSanityChecks(frame)) return frame ;
589
590 // Extract projection configuration from command list
591 RooCmdConfig pc("RooSimultaneous::plotOn(" + std::string(GetName()) + ")");
592 pc.defineString("sliceCatState","SliceCat",0,"",true) ;
593 pc.defineDouble("scaleFactor","Normalization",0,1.0) ;
594 pc.defineInt("scaleType","Normalization",0,RooAbsPdf::Relative) ;
595 pc.defineObject("sliceCatList","SliceCat",0,nullptr,true) ;
596 // This dummy is needed for plotOn to recognize the "SliceCatMany" command.
597 // It is not used directly, but the "SliceCat" commands are nested in it.
598 // Removing this dummy definition results in "ERROR: unrecognized command: SliceCatMany".
599 pc.defineObject("dummy1","SliceCatMany",0) ;
600 pc.defineSet("projSet","Project",0) ;
601 pc.defineSet("sliceSet","SliceVars",0) ;
602 pc.defineSet("projDataSet","ProjData",0) ;
603 pc.defineObject("projData","ProjData",1) ;
604 pc.defineMutex("Project","SliceVars") ;
605 pc.allowUndefined() ; // there may be commands we don't handle here
606
607 // Process and check varargs
608 pc.process(cmdList) ;
609 if (!pc.ok(true)) {
610 return frame ;
611 }
612
613 RooAbsData* projData = static_cast<RooAbsData*>(pc.getObject("projData")) ;
614 const RooArgSet* projDataSet = pc.getSet("projDataSet");
615 const RooArgSet* sliceSetTmp = pc.getSet("sliceSet") ;
616 std::unique_ptr<RooArgSet> sliceSet( sliceSetTmp ? (static_cast<RooArgSet*>(sliceSetTmp->Clone())) : nullptr );
617 const RooArgSet* projSet = pc.getSet("projSet") ;
618 double scaleFactor = pc.getDouble("scaleFactor") ;
619 ScaleType stype = (ScaleType) pc.getInt("scaleType") ;
620
621
622 // Look for category slice arguments and add them to the master slice list if found
623 const char* sliceCatState = pc.getString("sliceCatState",nullptr,true) ;
624 const RooLinkedList& sliceCatList = pc.getObjectList("sliceCatList") ;
625 if (sliceCatState) {
626
627 // Make the master slice set if it doesnt exist
628 if (!sliceSet) {
629 sliceSet = std::make_unique<RooArgSet>();
630 }
631
632 // Prepare comma separated label list for parsing
634
635 // Loop over all categories provided by (multiple) Slice() arguments
636 unsigned int tokenIndex = 0;
638 const char* slabel = tokenIndex >= catTokens.size() ? nullptr : catTokens[tokenIndex++].c_str();
639
640 if (slabel) {
641 // Set the slice position to the value indicated by slabel
642 scat->setLabel(slabel) ;
643 // Add the slice category to the master slice set
644 sliceSet->add(*scat,false) ;
645 }
646 }
647 }
648
649 // Check if we have a projection dataset
650 if (!projData) {
651 coutE(InputArguments) << "RooSimultaneous::plotOn(" << GetName() << ") ERROR: must have a projection dataset for index category" << std::endl ;
652 return frame ;
653 }
654
655 // Make list of variables to be projected
657 if (sliceSet) {
658 makeProjectionSet(frame->getPlotVar(),frame->getNormVars(),projectedVars,true) ;
659
660 // Take out the sliced variables
661 for (const auto sliceArg : *sliceSet) {
662 RooAbsArg* arg = projectedVars.find(sliceArg->GetName()) ;
663 if (arg) {
664 projectedVars.remove(*arg) ;
665 } else {
666 coutI(Plotting) << "RooAbsReal::plotOn(" << GetName() << ") slice variable "
667 << sliceArg->GetName() << " was not projected anyway" << std::endl ;
668 }
669 }
670 } else if (projSet) {
671 makeProjectionSet(frame->getPlotVar(),projSet,projectedVars,false) ;
672 } else {
673 makeProjectionSet(frame->getPlotVar(),frame->getNormVars(),projectedVars,true) ;
674 }
675
676 bool projIndex(false) ;
677
678 if (!_indexCat.arg().isDerived()) {
679 // *** Error checking for a fundamental index category ***
680 //cout << "RooSim::plotOn: index is fundamental" << std::endl ;
681
682 // Check that the provided projection dataset contains our index variable
683 if (!projData->get()->find(_indexCat.arg().GetName())) {
684 coutE(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") ERROR: Projection over index category "
685 << "requested, but projection data set doesn't contain index category" << std::endl ;
686 return frame ;
687 }
688
689 if (projectedVars.find(_indexCat.arg().GetName())) {
691 }
692
693 } else {
694 // *** Error checking for a composite index category ***
695
696 // Determine if any servers of the index category are in the projectedVars
698 bool anyServers(false) ;
699 for (const auto server : flattenedCatList()) {
700 if (projectedVars.find(server->GetName())) {
702 projIdxServers.add(*server) ;
703 }
704 }
705
706 // Check that the projection dataset contains all the
707 // index category components we're projecting over
708
709 // Determine if all projected servers of the index category are in the projection dataset
710 bool allServers(true) ;
711 std::string missing;
712 for (const auto server : projIdxServers) {
713 if (!projData->get()->find(server->GetName())) {
715 missing = server->GetName();
716 }
717 }
718
719 if (!allServers) {
720 coutE(Plotting) << "RooSimultaneous::plotOn(" << GetName()
721 << ") ERROR: Projection dataset doesn't contain complete set of index categories to do projection."
722 << "\n\tcategory " << missing << " is missing." << std::endl ;
723 return frame ;
724 }
725
726 if (anyServers) {
727 projIndex = true ;
728 }
729 }
730
731 // Calculate relative weight fractions of components
732 std::unique_ptr<Roo1DTable> wTable( projData->table(_indexCat.arg()) );
733
734 // Clone the index category to be able to cycle through the category states for plotting without
735 // affecting the category state of our instance
737 RooArgSet(*_indexCat).snapshot(idxCloneSet, true);
738 auto idxCatClone = static_cast<RooAbsCategoryLValue*>(idxCloneSet.find(_indexCat->GetName()) );
740
741 // Make list of category columns to exclude from projection data
742 std::unique_ptr<RooArgSet> idxCompSliceSet( idxCatClone->getObservables(frame->getNormVars()) );
743
744 // If we don't project over the index, just do the regular plotOn
745 if (!projIndex) {
746
747 coutI(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") plot on " << frame->getPlotVar()->GetName()
748 << " represents a slice in the index category (" << _indexCat.arg().GetName() << ")" << std::endl ;
749
750 // Reduce projData: take out fitCat (component) columns and entries that don't match selected slice
751 // Construct cut string to only select projection data event that match the current slice
752
753 // Make cut string to exclude rows from projection data
754 if (sliceSet) {
756 if (auto* slicedComponent = static_cast<const RooAbsCategory*>(sliceSet->find(*idxComp))) {
757 idxComp->setIndex(slicedComponent->getCurrentIndex(), false);
758 }
759 }
760 }
762
763 // Make temporary projData without RooSim index category components
764 RooArgSet projDataVars(*projData->get()) ;
765 projDataVars.remove(*idxCompSliceSet,true,true) ;
766
767 std::unique_ptr<RooAbsData>
769
770 // Override normalization and projection dataset
771 RooCmdArg tmp1 = RooFit::Normalization(scaleFactor*wTable->getFrac(idxCatClone->getCurrentLabel()),stype) ;
773
774 // WVE -- do not adjust normalization for asymmetry plots
776 if (!cmdList.find("Asymmetry")) {
777 cmdList2.Add(&tmp1) ;
778 }
779 cmdList2.Add(&tmp2) ;
780
781 // Plot single component
782 RooPlot* retFrame = getPdf(idxCatClone->getCurrentLabel())->plotOn(frame,cmdList2);
783 return retFrame ;
784 }
785
786 // If we project over the index, plot using a temporary RooAddPdf
787 // using the weights from the data as coefficients
788
789 // Build the list of indexCat components that are sliced
790 idxCompSliceSet->remove(projectedVars,true,true) ;
791
792 // Make a new expression that is the weighted sum of requested components
795//RooAbsPdf* pdf ;
796 double sumWeight(0) ;
798
799 idxCatClone->setLabel(proxy->name()) ;
800
801 // Determine if this component is the current slice (if we slice)
802 bool skip(false) ;
803 for (const auto idxSliceCompArg : *idxCompSliceSet) {
804 const auto idxSliceComp = static_cast<RooAbsCategory*>(idxSliceCompArg);
805 RooAbsCategory* idxComp = static_cast<RooAbsCategory*>(idxCloneSet.find(idxSliceComp->GetName())) ;
806 if (idxComp->getCurrentIndex()!=idxSliceComp->getCurrentIndex()) {
807 skip=true ;
808 break ;
809 }
810 }
811 if (skip) continue ;
812
813 // Instantiate a RRV holding this pdfs weight fraction
814 wgtCompList.addOwned(std::make_unique<RooRealVar>(proxy->name(),"coef",wTable->getFrac(proxy->name())));
815 sumWeight += wTable->getFrac(proxy->name()) ;
816
817 // Add the PDF to list list
818 pdfCompList.add(proxy->arg()) ;
819 }
820
822 RooAddPdf plotVar{plotVarName,"weighted sum of RS components",pdfCompList,wgtCompList};
823
824 // Fix appropriate coefficient normalization in plot function
825 if (!_plotCoefNormSet.empty()) {
826 plotVar.fixAddCoefNormalization(_plotCoefNormSet) ;
827 }
828
829 std::unique_ptr<RooAbsData> projDataTmp;
831 if (projData) {
832
833 // Construct cut string to only select projection data event that match the current slice
835
836 // Make temporary projData without RooSim index category components
837 RooArgSet projDataVars(*projData->get()) ;
839 _indexCat.arg().getObservables(frame->getNormVars(), idxCatServers) ;
840
841 projDataVars.remove(idxCatServers,true,true) ;
842
843 projDataTmp = std::unique_ptr<RooAbsData>{projData->reduce(RooFit::SelectVars(projDataVars), RooFit::Cut(cutString.c_str()))};
844
845
846
847 if (projSet) {
848 projSetTmp.add(*projSet) ;
849 projSetTmp.remove(idxCatServers,true,true);
850 }
851 }
852
853
854 if (_indexCat.arg().isDerived() && !idxCompSliceSet->empty()) {
855 coutI(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") plot on " << frame->getPlotVar()->GetName()
856 << " represents a slice in index category components " << *idxCompSliceSet << std::endl ;
857
859 _indexCat.arg().getObservables(frame->getNormVars(), idxCompProjSet) ;
860 idxCompProjSet.remove(*idxCompSliceSet,true,true) ;
861 if (!idxCompProjSet.empty()) {
862 coutI(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") plot on " << frame->getPlotVar()->GetName()
863 << " averages with data index category components " << idxCompProjSet << std::endl ;
864 }
865 } else {
866 coutI(Plotting) << "RooSimultaneous::plotOn(" << GetName() << ") plot on " << frame->getPlotVar()->GetName()
867 << " averages with data index category (" << _indexCat.arg().GetName() << ")" << std::endl ;
868 }
869
870
871 // Override normalization and projection dataset
873
874 RooCmdArg tmp1 = RooFit::Normalization(scaleFactor*sumWeight,stype) ;
876 // WVE -- do not adjust normalization for asymmetry plots
877 if (!cmdList.find("Asymmetry")) {
878 cmdList2.Add(&tmp1) ;
879 }
880 cmdList2.Add(&tmp2) ;
881
882 RooPlot* frame2 ;
883 if (!projSetTmp.empty()) {
884 // Plot temporary function
886 cmdList2.Add(&tmp3) ;
887 frame2 = plotVar.plotOn(frame,cmdList2) ;
888 } else {
889 // Plot temporary function
890 frame2 = plotVar.plotOn(frame,cmdList2) ;
891 }
892
893 return frame2 ;
894}
895
896
897////////////////////////////////////////////////////////////////////////////////
898/// Interface function used by test statistics to freeze choice of observables
899/// for interpretation of fraction coefficients. Needed here because a RooSimultaneous
900/// works like a RooAddPdf when plotted
901
907
908
909////////////////////////////////////////////////////////////////////////////////
910/// Interface function used by test statistics to freeze choice of range
911/// for interpretation of fraction coefficients. Needed here because a RooSimultaneous
912/// works like a RooAddPdf when plotted
913
918
919
920
921
922////////////////////////////////////////////////////////////////////////////////
923
925 const RooArgSet* auxProto, bool verbose, bool autoBinned, const char* binnedTag) const
926{
927 const char* idxCatName = _indexCat.arg().GetName() ;
928
929 if (vars.find(idxCatName) && prototype==nullptr
930 && (auxProto==nullptr || auxProto->empty())
931 && (autoBinned || (binnedTag && strlen(binnedTag)))) {
932
933 // Return special generator config that can also do binned generation for selected states
934 return new RooSimSplitGenContext(*this,vars,verbose,autoBinned,binnedTag) ;
935
936 } else {
937
938 // Return regular generator config ;
939 return genContext(vars,prototype,auxProto,verbose) ;
940 }
941}
942
943
944
945////////////////////////////////////////////////////////////////////////////////
946/// Return specialized generator context for simultaneous p.d.f.s
947
949 const RooArgSet* auxProto, bool verbose) const
950{
951 RooArgSet allVars{vars};
952 if(prototype) allVars.add(*prototype->get());
953
956
957 // Not generating index cat: we better error out because it's not clear what
958 // the user expects here. Does the user want to generate according to the
959 // currently-selected pdf? Or does the user want to generate global
960 // observable values according to the union of all category pdfs?
961 // Print an error and tell the user what to do to explicitly.
962 if(catsAmongAllVars.empty()) {
963 coutE(InputArguments) << "RooSimultaneous::generateSimGlobal(" << GetName()
964 << ") asking to generate without the index category!\n"
965 << "It's not clear what to do. you probably want to either:\n"
966 << "\n"
967 << " 1. Generate according to the currently-selected pdf.\n"
968 << " Please do this explicitly with:\n"
969 << " simpdf->getPdf(simpdf->indexCat().getCurrentLabel())->generate(vars, ...)\n"
970 << "\n"
971 << " 1. Generate global observable values according to the union of all component pdfs.\n"
972 << " For this, please use simpdf->generateSimGlobal(vars, ...)\n"
973 << std::endl;
974 return nullptr;
975 }
976
978 if(prototype) {
979 prototype->get()->selectCommon(flattenedCatList(), catsAmongProtoVars);
980
981 if(!catsAmongProtoVars.empty() && catsAmongProtoVars.size() != flattenedCatList().size()) {
982 // Abort if we have only part of the servers
983 coutE(Plotting) << "RooSimultaneous::genContext: ERROR: prototype must include either all "
984 << " components of the RooSimultaneous index category or none " << std::endl;
985 return nullptr;
986 }
987 }
988
989 return new RooSimGenContext(*this,vars,prototype,auxProto,verbose) ;
990}
991
992
993
994
995////////////////////////////////////////////////////////////////////////////////
996
998 const RooArgSet* nset,
999 double scaleFactor,
1001 bool showProgress) const
1002{
1003 if (RooAbsReal::fillDataHist (hist, nset, scaleFactor,
1004 correctForBinVolume, showProgress) == nullptr)
1005 return nullptr;
1006
1007 const double sum = hist->sumEntries();
1008 if (sum != 0) {
1009 for (int i=0 ; i<hist->numEntries() ; i++) {
1010 hist->set(i, hist->weight(i) / sum, 0.);
1011 }
1012 }
1013
1014 return hist;
1015}
1016
1017
1018
1019
1020////////////////////////////////////////////////////////////////////////////////
1021/// Special generator interface for generation of 'global observables' -- for RooStats tools.
1022///
1023/// \note Why one can't just use RooAbsPdf::generate()? That's because when
1024/// using the regular generate() method, a specific component pdf is selected
1025/// for each generated entry according to the index category value. However,
1026/// global observable values are independent of the current index category,
1027/// which can best be illustrated with the case where a global observable
1028/// corresponds to a nuisance parameter that is relevant for multiple channels.
1029/// So the interpretation of what is an entry in the generated dataset is very
1030/// different, hence the separate function.
1031
1033{
1034 // Generating the index category together with the global observables doesn't make any sense.
1037 if(!catsAmongAllVars.empty()) {
1038 coutE(InputArguments) << "RooSimultaneous::generateSimGlobal(" << GetName()
1039 << ") asking to generate global obserables at the same time as the index category!\n"
1040 << "This doesn't make any sense: global observables are generally not related to a specific channel.\n"
1041 << std::endl;
1042 return nullptr;
1043 }
1044
1045 // Make set with clone of variables (placeholder for output)
1047 whatVars.snapshot(globClone);
1048
1049 auto data = std::make_unique<RooDataSet>("gensimglobal","gensimglobal",whatVars);
1050
1051 for (Int_t i=0 ; i<nEvents ; i++) {
1052 for (const auto& nameIdx : indexCat()) {
1053
1054 // Get pdf associated with state from simpdf
1055 RooAbsPdf* pdftmp = getPdf(nameIdx.first);
1056
1058 pdftmp->getObservables(&whatVars, globtmp) ;
1059
1060 // If there are any, generate only global variables defined by the pdf
1061 // associated with this state and transfer values to output placeholder.
1062 if (!globtmp.empty()) {
1063 globClone.assign(*std::unique_ptr<RooDataSet>{pdftmp->generate(globtmp,1)}->get(0)) ;
1064 }
1065 }
1066 data->add(globClone) ;
1067 }
1068
1069 return RooFit::makeOwningPtr(std::move(data));
1070}
1071
1072
1073/// Wraps the components of this RooSimultaneous in RooBinSamplingPdfs.
1074/// \param[in] data The dataset to be used in the eventual fit, used to figure
1075/// out the observables and whether the dataset is binned.
1076/// \param[in] precision Precision argument for all created RooBinSamplingPdfs.
1078
1079 if (precision < 0.) return;
1080
1082
1083 for (auto const &item : this->indexCat()) {
1084
1085 auto const &catName = item.first;
1086 auto &pdf = *this->getPdf(catName);
1087
1088 if (auto newSamplingPdf = RooBinSamplingPdf::create(pdf, data, precision)) {
1089 // Set the "ORIGNAME" attribute the indicate to
1090 // RooAbsArg::redirectServers() which pdf should be replaced by this
1091 // RooBinSamplingPdf in the RooSimultaneous.
1092 newSamplingPdf->setAttribute(
1093 (std::string("ORIGNAME:") + pdf.GetName()).c_str());
1094 newSamplingPdfs.addOwned(std::move(newSamplingPdf));
1095 }
1096 }
1097
1098 this->redirectServers(newSamplingPdfs, false, true);
1099 this->addOwnedComponents(std::move(newSamplingPdfs));
1100}
1101
1102
1103/// Wraps the components of this RooSimultaneous in RooBinSamplingPdfs, with a
1104/// different precision parameter for each component.
1105/// \param[in] data The dataset to be used in the eventual fit, used to figure
1106/// out the observables and whether the dataset is binned.
1107/// \param[in] precisions The map that gives the precision argument for each
1108/// component in the RooSimultaneous. The keys are the pdf names. If
1109/// there is no value for a given component, it will not use the bin
1110/// integration. Otherwise, the value has the same meaning than in
1111/// the IntegrateBins() command argument for RooAbsPdf::fitTo().
1112/// \param[in] useCategoryNames If this flag is set, the category names will be
1113/// used to look up the precision in the precisions map instead of
1114/// the pdf names.
1116 std::map<std::string, double> const& precisions,
1117 bool useCategoryNames /*=false*/) {
1118
1119 constexpr double defaultPrecision = -1.;
1120
1122
1123 for (auto const &item : this->indexCat()) {
1124
1125 auto const &catName = item.first;
1126 auto &pdf = *this->getPdf(catName);
1127 std::string pdfName = pdf.GetName();
1128
1129 auto found = precisions.find(useCategoryNames ? catName : pdfName);
1130 const double precision =
1131 found != precisions.end() ? found->second : defaultPrecision;
1132 if (precision < 0.)
1133 continue;
1134
1135 if (auto newSamplingPdf = RooBinSamplingPdf::create(pdf, data, precision)) {
1136 // Set the "ORIGNAME" attribute the indicate to
1137 // RooAbsArg::redirectServers() which pdf should be replaced by this
1138 // RooBinSamplingPdf in the RooSimultaneous.
1139 newSamplingPdf->setAttribute(
1140 (std::string("ORIGNAME:") + pdf.GetName()).c_str());
1141 newSamplingPdfs.addOwned(std::move(newSamplingPdf));
1142 }
1143 }
1144
1145 this->redirectServers(newSamplingPdfs, false, true);
1146 this->addOwnedComponents(std::move(newSamplingPdfs));
1147}
1148
1149/// Internal utility function to get a list of all category components for this
1150/// RooSimultaneous. The output contains only the index category if it is a
1151/// RooCategory, or the list of all category components if it is a
1152/// RooSuperCategory.
1154{
1155 // Note that the index category of a RooSimultaneous can only be of type
1156 // RooCategory or RooSuperCategory, because these are the only classes that
1157 // inherit from RooAbsCategoryLValue.
1158 if (auto superCat = dynamic_cast<RooSuperCategory const*>(&_indexCat.arg())) {
1159 return superCat->inputCatList();
1160 }
1161
1162 if(!_indexCatSet) {
1163 _indexCatSet = std::make_unique<RooArgSet>(_indexCat.arg());
1164 }
1165 return *_indexCatSet;
1166}
1167
1168namespace {
1169
1170void markObs(RooAbsArg *arg, std::string const &prefix, RooArgSet const &normSet)
1171{
1172 for (RooAbsArg *server : arg->servers()) {
1173 if (server->isFundamental() && normSet.find(*server)) {
1174 markObs(server, prefix, normSet);
1175 server->setAttribute("__obs__");
1176 } else if (!server->isFundamental()) {
1177 markObs(server, prefix, normSet);
1178 }
1179 }
1180}
1181
1182void prefixArgs(RooAbsArg *arg, std::string const &prefix, RooArgSet const &normSet)
1183{
1184 if (!arg->getStringAttribute("__prefix__")) {
1185 arg->SetName((prefix + arg->GetName()).c_str());
1186 arg->setStringAttribute("__prefix__", prefix.c_str());
1187 }
1188 for (RooAbsArg *server : arg->servers()) {
1189 if (server->isFundamental() && normSet.find(*server)) {
1190 prefixArgs(server, prefix, normSet);
1191 } else if (!server->isFundamental()) {
1192 prefixArgs(server, prefix, normSet);
1193 }
1194 }
1195}
1196
1197} // namespace
1198
1199std::unique_ptr<RooAbsArg>
1201{
1202 std::unique_ptr<RooSimultaneous> newSimPdf{static_cast<RooSimultaneous *>(this->Clone())};
1203
1204 const char *rangeName = this->getStringAttribute("RangeName");
1205 bool splitRange = this->getAttribute("SplitRange");
1206
1208 std::vector<std::string> catNames;
1209
1210 for (auto *proxy : static_range_cast<RooRealProxy *>(newSimPdf->_pdfProxyList)) {
1211 catNames.emplace_back(proxy->GetName());
1212 std::string const &catName = catNames.back();
1213 const std::string prefix = "_" + catName + "_";
1214
1215 const std::string origname = proxy->arg().GetName();
1216
1217 auto pdfClone = RooHelpers::cloneTreeWithSameParameters(static_cast<RooAbsPdf const &>(proxy->arg()), &normSet);
1218
1219 markObs(pdfClone.get(), prefix, normSet);
1220
1221 std::unique_ptr<RooArgSet> pdfNormSet{
1222 std::unique_ptr<RooArgSet>(pdfClone->getVariables())->selectByAttrib("__obs__", true)};
1223 std::unique_ptr<RooArgSet> condVarSet{
1224 std::unique_ptr<RooArgSet>(pdfClone->getVariables())->selectByAttrib("__conditional__", true)};
1225
1226 pdfNormSet->remove(*condVarSet, true, true);
1227
1228 if (rangeName) {
1230 }
1231
1233 pdfContext.setLikelihoodMode(ctx.likelihoodMode());
1234 auto *pdfFinal = pdfContext.compile(*pdfClone, *newSimPdf, *pdfNormSet);
1235
1236 // We can only prefix the observables after everything related the
1237 // compiling of the compute graph for the normalization set is done. This
1238 // is because of a subtlety in conditional RooProdPdfs, which stores the
1239 // normalization sets for the individual pdfs in RooArgSets that are
1240 // disconnected from the computation graph, so we have no control over
1241 // them. An alternative would be to use recursive server re-direction,
1242 // but this has more performance overhead.
1243 prefixArgs(pdfFinal, prefix, normSet);
1244
1245 pdfFinal->fixAddCoefNormalization(*pdfNormSet, false);
1246
1247 pdfClone->SetName((std::string("_") + pdfClone->GetName()).c_str());
1248 pdfFinal->addOwnedComponents(std::move(pdfClone));
1249
1250 pdfFinal->setAttribute(("ORIGNAME:" + origname).c_str());
1251 newPdfs.add(*pdfFinal);
1252
1253 // We will remove the old pdf server because we will fill the new ones by
1254 // hand via the creation of new proxies.
1255 newSimPdf->removeServer(const_cast<RooAbsReal &>(proxy->arg()), true);
1256 }
1257
1258 // Replace pdfs with compiled pdfs. Don't use RooAbsArg::redirectServers()
1259 // here, because it doesn't support replacing two servers with the same name
1260 // (it can happen in a RooSimultaneous that two pdfs have the same name).
1261
1262 // First delete old proxies (we have already removed the servers before).
1263 newSimPdf->_pdfProxyList.Delete();
1264
1265 // Recreate the _pdfProxyList with the compiled pdfs
1266 for (std::size_t i = 0; i < newPdfs.size(); ++i) {
1267 const char *label = catNames[i].c_str();
1268 newSimPdf->_pdfProxyList.Add(
1269 new RooRealProxy(label, label, newSimPdf.get(), *static_cast<RooAbsReal *>(newPdfs[i])));
1270 }
1271
1272 ctx.compileServers(*newSimPdf, normSet); // to trigger compiling also the index category
1273
1274 return newSimPdf;
1275}
size_t size(const MatrixT &matrix)
retrieve the size of a square matrix
#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.
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
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
const_iterator begin() const
const_iterator end() const
Common abstract base class for objects that represent a value and a "shape" in RooFit.
Definition RooAbsArg.h:76
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.
void setStringAttribute(const Text_t *key, const Text_t *value)
Associate string 'value' to this object under key 'key'.
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.
bool getAttribute(const Text_t *name) const
Check if a named attribute is set. By default, all attributes are unset.
TObject * Clone(const char *newname=nullptr) const override
Make a clone of an object using the Streamer facility.
Definition RooAbsArg.h:88
Abstract base class for objects that represent a discrete value that can be set from the outside,...
A space to attach TBranches.
bool contains(const RooAbsArg &var) const
Check if collection contains an argument with the same name as var.
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
virtual bool addOwned(RooAbsArg &var, bool silent=false)
Add an argument and transfer the ownership to the collection.
RooAbsArg * find(const char *name) const
Find object with given name in list.
Abstract base class for binned and unbinned datasets.
Definition RooAbsData.h:57
virtual const RooArgSet * get() const
Definition RooAbsData.h:101
virtual Roo1DTable * table(const RooArgSet &catSet, const char *cuts="", const char *opts="") const
Construct table for product of categories in catSet.
RooFit::OwningPtr< RooAbsData > reduce(const RooCmdArg &arg1, const RooCmdArg &arg2={}, const RooCmdArg &arg3={}, const RooCmdArg &arg4={}, const RooCmdArg &arg5={}, const RooCmdArg &arg6={}, const RooCmdArg &arg7={}, const RooCmdArg &arg8={}) const
Create a reduced copy of this dataset.
virtual Int_t numEntries() const
Return number of entries in dataset, i.e., count unweighted entries.
Abstract base class for generator contexts of RooAbsPdf objects.
Abstract interface for all probability density functions.
Definition RooAbsPdf.h:32
RooArgSet const * _normSet
Normalization integral (owned by _normMgr)
Definition RooAbsPdf.h:316
RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1={}, const RooCmdArg &arg2={}, const RooCmdArg &arg3={}, const RooCmdArg &arg4={}, const RooCmdArg &arg5={}, const RooCmdArg &arg6={}, const RooCmdArg &arg7={}, const RooCmdArg &arg8={}, const RooCmdArg &arg9={}, const RooCmdArg &arg10={}) const override
Helper calling plotOn(RooPlot*, RooLinkedList&) const.
Definition RooAbsPdf.h:116
bool canBeExtended() const
If true, PDF can provide extended likelihood term.
Definition RooAbsPdf.h:214
@ CanBeExtended
Definition RooAbsPdf.h:208
@ MustBeExtended
Definition RooAbsPdf.h:208
@ CanNotBeExtended
Definition RooAbsPdf.h:208
Abstract base class for objects that represent a real value and implements functionality common to al...
Definition RooAbsReal.h:63
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:107
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...
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:24
RooArgSet * snapshot(bool deepCopy=true) const
Use RooAbsCollection::snapshot(), but return as RooArgSet.
Definition RooArgSet.h:159
RooArgSet * selectCommon(const RooAbsCollection &refColl) const
Use RooAbsCollection::selecCommon(), but return as RooArgSet.
Definition RooArgSet.h:154
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.
Named container for two doubles, two integers two object points and three string pointers that can be...
Definition RooCmdArg.h:26
Configurable parser for RooCmdArg named arguments.
void defineMutex(const char *head, Args_t &&... tail)
Define arguments where any pair is mutually exclusive.
bool process(const RooCmdArg &arg)
Process given RooCmdArg.
double getDouble(const char *name, double defaultValue=0.0) const
Return double property registered with name 'name'.
bool defineDouble(const char *name, const char *argName, int doubleNum, double defValue=0.0)
Define double property name 'name' mapped to double in slot 'doubleNum' in RooCmdArg with name argNam...
RooArgSet * getSet(const char *name, RooArgSet *set=nullptr) const
Return RooArgSet property registered with name 'name'.
bool defineSet(const char *name, const char *argName, int setNum, const RooArgSet *set=nullptr)
Define TObject property name 'name' mapped to object in slot 'setNum' in RooCmdArg with name argName ...
bool ok(bool verbose) const
Return true of parsing was successful.
bool defineObject(const char *name, const char *argName, int setNum, const TObject *obj=nullptr, bool isArray=false)
Define TObject property name 'name' mapped to object in slot 'setNum' in RooCmdArg with name argName ...
const char * getString(const char *name, const char *defaultValue="", bool convEmptyToNull=false) const
Return string property registered with name 'name'.
bool defineString(const char *name, const char *argName, int stringNum, const char *defValue="", bool appendMode=false)
Define double property name 'name' mapped to double in slot 'stringNum' in RooCmdArg with name argNam...
const RooLinkedList & getObjectList(const char *name) const
Return list of objects registered with name 'name'.
bool defineInt(const char *name, const char *argName, int intNum, int defValue=0)
Define integer property name 'name' mapped to integer in slot 'intNum' in RooCmdArg with name argName...
void allowUndefined(bool flag=true)
If flag is true the processing of unrecognized RooCmdArgs is not considered an error.
int getInt(const char *name, int defaultValue=0) const
Return integer property registered with name 'name'.
TObject * getObject(const char *name, TObject *obj=nullptr) const
Return TObject property registered with name 'name'.
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...
Container class to hold N-dimensional binned data.
Definition RooDataHist.h:40
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.
Container class to hold unbinned data.
Definition RooDataSet.h:32
void compileServers(RooAbsArg &arg, RooArgSet const &normSet)
Collection class for internal use, storing a collection of RooAbsArg pointers in a doubly linked list...
static const TNamed * ptr(const char *stringPtr)
Return a unique TNamed pointer for given C++ string.
Plot frame and a container for graphics objects within that frame.
Definition RooPlot.h:43
const RooArgSet * getNormVars() const
Definition RooPlot.h:146
RooAbsRealLValue * getPlotVar() const
Definition RooPlot.h:137
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
~RooSimultaneous() override
Destructor.
RooFit::OwningPtr< RooDataSet > generateSimGlobal(const RooArgSet &whatVars, Int_t nEvents) override
Special generator interface for generation of 'global observables' – for RooStats tools.
RooArgSet const & flattenedCatList() const
Internal utility function to get a list of all category components for this RooSimultaneous.
RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1, const RooCmdArg &arg2={}, const RooCmdArg &arg3={}, const RooCmdArg &arg4={}, const RooCmdArg &arg5={}, const RooCmdArg &arg6={}, const RooCmdArg &arg7={}, const RooCmdArg &arg8={}, const RooCmdArg &arg9={}, const RooCmdArg &arg10={}) const override
ExtendMode extendMode() const override
Examine the pdf components and check if one of them can be extended or must be extended.
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-style string that can also be constructed from a std::strin...
Joins several RooAbsCategoryLValue objects into a single category.
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:576
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:468
const char * GetName() const override
Returns name of object.
Definition TNamed.h:49
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:138
RooCmdArg SelectVars(const RooArgSet &vars)
RooCmdArg ProjWData(const RooAbsData &projData, bool binData=false)
RooCmdArg Project(const RooArgSet &projSet)
RooCmdArg Normalization(double scaleFactor)
RooCmdArg Cut(const char *cutSpec)
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::string makeSliceCutString(RooArgSet const &sliceDataSet)
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Definition CodegenImpl.h:67
T * OwningPtr
An alias for raw pointers for indicating that the return type of a RooFit function is an owning point...
Definition Config.h:35
OwningPtr< T > makeOwningPtr(std::unique_ptr< T > &&ptr)
Internal helper to turn a std::unique_ptr<T> into an OwningPtr.
Definition Config.h:40
std::unique_ptr< T > cloneTreeWithSameParameters(T const &arg, RooArgSet const *observables=nullptr)
Clone RooAbsArg object and reattach to original parameters.
std::string getRangeNameForSimComponent(std::string const &rangeName, bool splitRange, std::string const &catName)
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:2340