Logo ROOT  
Reference Guide
RooSimPdfBuilder.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///
19/// \class RooSimPdfBuilder
20///
21/// <b>This tool has now been superseded by RooSimWSTool</b>
22///
23/// <p>
24/// <tt>RooSimPdfBuilder</tt> is a powerful tool to build <tt>RooSimultaneous</tt>
25/// PDFs that are defined in terms component PDFs that are identical in
26/// structure, but have different parameters.
27/// </p>
28///
29/// <h2>Example</h2>
30///
31/// <p>
32/// The following example demonstrates the essence of <tt>RooSimPdfBuilder</tt>:
33/// Given a dataset D with a <tt>RooRealVar X</tt> and a <tt>RooCategory C</tt> that has
34/// state C1 and C2.
35/// <ul>
36/// <li> We want to fit the distribution of <tt>X</tt> with a Gaussian+ArgusBG PDF,
37/// <li> We want to fit the data subsets <tt>D(C==C1)</tt> and <tt>D(C==C2)</tt> separately and simultaneously.
38/// <li> The PDFs to fit data subsets D_C1 and D_C2 are identical except for
39/// <ul>
40/// <li> the kappa parameter of the ArgusBG PDF and
41/// <li> the sigma parameter of the gaussian PDF
42/// </ul>
43/// where each PDF will have its own copy of the parameter
44/// </ul>
45/// </p>
46/// <p>
47/// Coding this example directly with RooFit classes gives
48/// (we assume dataset D and variables C and X have been declared previously)
49/// </p>
50/// <pre>
51/// RooRealVar m("m","mean of gaussian",-10,10) ;
52/// RooRealVar s_C1("s_C1","sigma of gaussian C1",0,20) ;
53/// RooRealVar s_C2("s_C2","sigma of gaussian C2",0,20) ;
54/// RooGaussian gauss_C1("gauss_C1","gaussian C1",X,m,s_C1) ;
55/// RooGaussian gauss_C2("gauss_C2","gaussian C2",X,m,s_C2) ;
56///
57/// RooRealVar k_C1("k_C1","ArgusBG kappa parameter C1",-50,0) ;
58/// RooRealVar k_C2("k_C2","ArgusBG kappa parameter C2",-50,0) ;
59/// RooRealVar xm("xm","ArgusBG cutoff point",5.29) ;
60/// RooArgusBG argus_C1("argus_C1","argus background C1",X,k_C1,xm) ;
61/// RooArgusBG argus_C2("argus_C2","argus background C2",X,k_C2,xm) ;
62///
63/// RooRealVar gfrac("gfrac","fraction of gaussian",0.,1.) ;
64/// RooAddPdf pdf_C1("pdf_C1","gauss+argus_C1",RooArgList(gauss_C1,argus_C1),gfrac) ;
65/// RooAddPdf pdf_C2("pdf_C2","gauss+argus_C2",RooArgList(gauss_C2,argus_C2),gfrac) ;
66///
67/// RooSimultaneous simPdf("simPdf","simPdf",C) ;
68/// simPdf.addPdf(pdf_C1,"C1") ;
69/// simPdf.addPdf(pdf_C2,"C2") ;
70/// </pre>
71/// <p>
72/// Coding this example with RooSimPdfBuilder gives
73/// </p>
74/// <pre>
75/// RooRealVar m("m","mean of gaussian",-10,10) ;
76/// RooRealVar s("s","sigma of gaussian",0,20) ;
77/// RooGaussian gauss("gauss","gaussian",X,m,s) ;
78///
79/// RooRealVar k("k","ArgusBG kappa parameter",-50,0) ;
80/// RooRealVar xm("xm","ArgusBG cutoff point",5.29) ;
81/// RooArgusBG argus("argus","argus background",X,k,xm) ;
82///
83/// RooRealVar gfrac("gfrac","fraction of gaussian",0.,1.) ;
84/// RooAddPdf pdf("pdf","gauss+argus",RooArgList(gauss,argus),gfrac) ;
85///
86/// RooSimPdfBuilder builder(pdf) ;
87/// RooArgSet* cfg = builder.createProtoBuildConfig() ;
88/// dynamic_cast<RooStringVar&>((*cfg)["physModels"]) = "pdf" ; /// Name of the PDF we are going to work with
89/// dynamic_cast<RooStringVar&>((*cfg)["splitCats"]) = "C" ; /// Category used to differentiate sub-datasets
90/// dynamic_cast<RooStringVar&>((*cfg)["pdf"]) = "C : k,s" ; /// Prescription to taylor PDF parameters k and s
91/// /// for each data subset designated by C states
92/// RooSimultaneous* simPdf = builder.buildPdf(*cfg,&D) ;
93/// </pre>
94/// <p>
95/// The above snippet of code demonstrates the concept of <tt>RooSimPdfBuilder</tt>:
96/// the user defines a single <i>'prototype' PDF</i> that defines the structure of all
97/// PDF components of the <tt>RooSimultaneous</tt> PDF to be built. <tt>RooSimPdfBuilder</tt>
98/// then takes this prototype and replicates it as a component
99/// PDF for each state of the C index category.
100/// </p>
101/// <p>
102/// In the above example <tt>RooSimPdfBuilder</tt>
103/// will first replicate <tt>k</tt> and <tt>s</tt> into
104/// <tt>k_C1,k_C2</tt> and <tt>s_C1,s_C2</tt>, as prescribed in the
105/// configuration. Then it will recursively replicate all PDF nodes that depend on
106/// the 'split' parameter nodes: <tt>gauss</tt> into <tt>gauss_C1,C2</tt>, <tt>argus</tt>
107/// into <tt>argus_C1,C2</tt> and finally <tt>pdf</tt> into <tt>pdf_C1,pdf_C2</tt>.
108/// When PDFs for all states of C have been replicated
109/// they are assembled into a <tt>RooSimultaneous</tt> PDF, which is returned by the <tt>buildPdf()</tt>
110/// method.
111/// </p>
112/// <p>
113/// Although in this very simple example the use of <tt>RooSimPdfBuilder</tt> doesn't
114/// reduce the amount of code much, it is already easier to read and maintain
115/// because there is no duplicate code. As the complexity of the <tt>RooSimultaneous</tt>
116/// to be built increases, the advantages of <tt>RooSimPdfBuilder</tt> will become more and
117/// more apparent.
118/// </p>
119///
120///
121/// <h2>Builder configuration rules for a single prototype PDF</h2>
122/// <p>
123/// Each builder configuration needs at minumum two lines, <tt>physModels</tt> and <tt>splitCats</tt>, which identify
124/// the ingredients of the build. In this section we only explain the building rules for
125/// builds from a single prototype PDF. In that case the <tt>physModels</tt> line always reads
126/// </p>
127/// <pre>
128/// physModels = {pdfName}
129/// </pre>
130/// <p>
131/// The second line, <tt>splitCats</tt>, indicates which categories are going to be used to
132/// differentiate the various subsets of the 'master' input data set. You can enter
133/// a single category here, or multiple if necessary:
134/// </p>
135/// <pre>
136/// splitCats = {catName} [{catName} ...]
137/// </pre>
138/// <p>
139/// All listed splitcats must be <tt>RooCategories</tt> that appear in the dataset provided to
140/// <tt>RooSimPdfBuilder::buildPdf()</tt>
141/// </p>
142/// <p>
143/// The parameter splitting prescriptions, the essence of each build configuration
144/// can be supplied in a third line carrying the name of the pdf listed in <tt>physModels</tt>
145/// </p>
146/// <pre>
147/// pdfName = {splitCat} : {parameter} [,{parameter},....]
148/// </pre>
149/// <p>
150/// Each pdf can have only one line with splitting rules, but multiple rules can be
151/// supplied in each line, e.g.
152/// </p>
153/// <pre>
154/// pdfName = {splitCat} : {parameter} [,{parameter},....]
155/// {splitCat} : {parameter} [,{parameter},....]
156/// </pre>
157/// <p>
158/// Conversely, each parameter can only have one splitting prescription, but it may be split
159/// by multiple categories, e.g.
160/// </p>
161/// <pre>
162/// pdfName = {splitCat1},{splitCat2} : {parameter}
163/// </pre>
164/// <p>
165/// instructs <tt>RooSimPdfBuilder</tt> to build a <tt>RooSuperCategory</tt>
166/// of <tt>{splitCat1}</tt> and <tt>{splitCat2}</tt>
167/// and split <tt>{parameter}</tt> with that <tt>RooSuperCategory</tt>
168/// </p>
169/// <p>
170/// Here is an example of a builder configuration that uses several of the options discussed
171/// above:
172/// </p>
173/// <pre>
174/// physModels = pdf
175/// splitCats = tagCat runBlock
176/// pdf = tagCat : signalRes,bkgRes
177/// runBlock : fudgeFactor
178/// tagCat,runBlock : kludgeParam
179/// </pre>
180///
181/// <h2>How to enter configuration data</h2>
182///
183/// <p>
184/// The prototype builder configuration returned by
185/// <tt>RooSimPdfBuilder::createProtoBuildConfig()</tt> is a pointer to a <tt>RooArgSet</tt> filled with
186/// initially blank <tt>RooStringVars</tt> named <tt>physModels,splitCats</tt> and one additional for each
187/// PDF supplied to the <tt>RooSimPdfBuilders</tt> constructor (with the same name)
188/// </p>
189/// <p>
190/// In macro code, the easiest way to assign new values to these <tt>RooStringVars</tt>
191/// is to use <tt>RooArgSet</tt>s array operator and the <tt>RooStringVar</tt>s assignment operator, e.g.
192/// </p>
193/// <pre>
194/// (*config)["physModels"] = "Blah" ;
195/// </pre>
196/// <p>
197/// To enter multiple splitting rules simply separate consecutive rules by whitespace
198/// (not newlines), e.g.
199/// </p>
200/// <pre>
201/// (*config)["physModels"] = "Blah " /// << note trailing space here
202/// "Blah 2" ;
203/// </pre>
204/// <p>
205/// In this example, the C++ compiler will concatenate the two string literals (without inserting
206/// any whitespace), so the extra space after 'Blah' is important here.
207/// </p>
208/// <p>
209/// Alternatively, you can read the configuration from an ASCII file, as you can
210/// for any <tt>RooArgSet</tt> using <tt>RooArgSet::readFromFile()</tt>. In that case the ASCII file
211/// can follow the syntax of the examples above and the '<tt>\</tt>' line continuation
212/// sequence can be used to fold a long splitting rule over multiple lines.
213/// </p>
214/// <pre>
215/// RooArgSet* config = builder.createProtoBuildConfig() ;
216/// config->readFromFile("config.txt") ;
217///
218/// --- config.txt ----------------
219/// physModels = pdf
220/// splitCats = tagCat
221/// pdf = tagCat : bogusPar
222/// -------------------------------
223/// </pre>
224///
225///
226/// <h2>Working with multiple prototype PDFs</h2>
227/// <p>
228/// It is also possible to build a <tt>RooSimultaneous</tt> PDF from multiple PDF prototypes.
229/// This is appropriate for cases where the input prototype PDF would otherwise be
230/// a <tt>RooSimultaneous</tt> PDF by itself. In such cases we don't feed a single
231/// <tt>RooSimultaneous</tt> PDF into <tt>RooSimPdfBuilder</tt>, instead we feed it its ingredients and
232/// add a prescription to the builder configuration that corresponds to the
233/// PDF-category state mapping of the prototype <tt>RooSimultaneous</tt>.
234/// </p>
235/// <p>
236/// The constructor of the <tt>RooSimPdfBuilder</tt> will look as follows:
237/// </p>
238/// <pre>
239/// RooSimPdfBuilder builder(RooArgSet(pdfA,pdfB,...)) ;
240/// </pre>
241/// <p>
242/// The <tt>physModels</tt> line is now expanded to carry the pdf->state mapping information
243/// that the prototype <tt>RooSimultaneous</tt> would have. I.e.
244/// </p>
245/// <pre>
246/// physModels = mode : pdfA=modeA pdfB=modeB
247/// </pre>
248/// <p>
249/// is equivalent to a prototype <tt>RooSimultaneous</tt> constructed as
250/// </p>
251/// <pre>
252/// RooSimultanous simPdf("simPdf","simPdf",mode);
253/// simPdf.addPdf(pdfA,"modeA") ;
254/// simPdf.addPdf(pdfB,"modeB") ;
255/// </pre>
256/// <p>
257/// The rest of the builder configuration works the same, except that
258/// each prototype PDF now has its own set of splitting rules, e.g.
259/// </p>
260/// <pre>
261/// physModels = mode : pdfA=modeA pdfB=modeB
262/// splitCats = tagCat
263/// pdfA = tagCat : bogusPar
264/// pdfB = tagCat : fudgeFactor
265/// </pre>
266/// <p>
267/// Please note that
268/// <ul>
269/// <li> The master index category ('mode' above) doesn't have to be listed in
270/// <tt>splitCats</tt>, this is implicit.
271///
272/// <li> The number of splitting prescriptions goes by the
273/// number of prototype PDFs and not by the number of states of the
274/// master index category (mode in the above and below example).
275/// </ul>
276///
277/// In the following case:
278///</p>
279/// <pre>
280/// physModels = mode : pdfA=modeA pdfB=modeB pdfA=modeC pdfB=modeD
281/// </pre>
282/// <p>
283/// there are still only 2 sets of splitting rules: one for <tt>pdfA</tt> and one
284/// for <tt>pdfB</tt>. However, you <i>can</i> differentiate between <tt>modeA</tt> and <tt>modeC</tt> in
285/// the above example. The technique is to use <tt>mode</tt> as splitting category, e.g.
286/// </p>
287/// <pre>
288/// physModels = mode : pdfA=modeA pdfB=modeB pdfA=modeC pdfB=modeD
289/// splitCats = tagCat
290/// pdfA = tagCat : bogusPar
291/// mode : funnyPar
292/// pdfB = mode : kludgeFactor
293/// </pre>
294/// <p>
295/// will result in an individual set of <tt>funnyPar</tt> parameters for <tt>modeA</tt> and <tt>modeC</tt>
296/// labeled <tt>funnyPar_modeA</tt> and <tt>funnyPar_modeB</tt> and an individual set of
297/// kludgeFactor parameters for <tt>pdfB</tt>, <tt>kludgeFactor_modeB</tt> and <tt>kludgeFactor_modeD</tt>.
298/// Please note that for splits in the master index category (mode) only the
299/// applicable states are built (A,C for <tt>pdfA</tt>, B,D for <tt>pdfB</tt>)
300/// </p>
301///
302///
303/// <h2>Advanced options</h2>
304///
305/// <h4>Partial splits</h4>
306/// <p>
307/// You can request to limit the list of states of each splitCat that
308/// will be considered in the build. This limitation is requested in the
309/// each build as follows:
310/// </p>
311/// <pre>
312/// splitCats = tagCat(Lep,Kao) RunBlock(Run1)
313/// </pre>
314/// <p>
315/// In this example the splitting of <tt>tagCat</tt> is limited to states <tt>Lep,Kao</tt>
316/// and the splitting of <tt>runBlock</tt> is limited to <tt>Run1</tt>. The splits apply
317/// globally to each build, i.e. every parameter split requested in this
318/// build will be limited according to these specifications.
319/// </p>
320/// <p>
321/// NB: Partial builds have no pdf associated with the unbuilt states of the
322/// limited splits. Running such a pdf on a dataset that contains data with
323/// unbuilt states will result in this data being ignored completely.
324/// </p>
325///
326///
327/// <h4>Non-trivial splits</h4>
328/// <p>
329/// It is possible to make non-trivial parameter splits with <tt>RooSimPdfBuilder</tt>.
330/// Trivial splits are considered simple splits in one (fundamental) category
331/// in the dataset or a split in a <tt>RooSuperCategory</tt> 'product' of multiple
332/// fundamental categories in the dataset. Non-trivial splits can be performed
333/// using an intermediate 'category function' (<tt>RooMappedCategory,
334/// RooGenericCategory,RooThresholdCategory</tt> etc), i.e. any <tt>RooAbsCategory</tt>
335/// derived objects that calculates its output as function of one or more
336/// input <tt>RooRealVars</tt> and/or <tt>RooCategories</tt>.
337/// </p>
338/// <p>
339/// Such 'function categories' objects must be constructed by the user prior
340/// to building the PDF. In the <tt>RooSimPdfBuilder::buildPdf()</tt> function these
341/// objects can be passed in an optional <tt>RooArgSet</tt> called 'auxiliary categories':
342/// </p>
343/// <pre>
344/// const <tt>RooSimultaneous</tt>* buildPdf(const RooArgSet& buildConfig, const RooAbsData* dataSet,
345/// const RooArgSet& auxSplitCats, bool verbose=false) {
346/// ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
347/// </pre>
348/// <p>
349/// Objects passed in this argset can subsequently be used in the build configuration, e.g.
350/// </p>
351/// <pre>
352/// RooMappedCategory tagMap("tagMap","Mapped tagging category",tagCat,"CutBased") ;
353/// tagMap.map("Lep","CutBased") ;
354/// tagMap.map("Kao","CutBased") ;
355/// tagMap.map("NT*","NeuralNet") ;
356/// ...
357/// builder.buildPdf(config,D,tagMap) ;
358/// ^^^^^^
359///<Contents of config>
360/// physModels = pdf
361/// splitCats = tagCat runBlock
362/// pdf = tagCat : signalRes
363/// tagMap : fudgeFactor
364/// ^^^^^^
365/// </pre>
366/// <p>
367/// In the above example <tt>signalRes</tt> will be split in <tt>signalRes_Kao,signalRes_Lep,
368/// signalRes_NT1,signalRes_NT2</tt>, while <tt>fudgeFactor</tt> will be split in <tt>fudgeFactor_CutBased</tt>
369/// and <tt>fudgeFactor_NeuralNet</tt>.
370/// </p>
371/// <p>
372/// Category functions passed in the auxSplitCats <tt>RooArgSet</tt> can be used regularly
373/// in the splitting configuration. They should not be listed in <tt>splitCats</tt>,
374/// but must be able to be expressed <i>completely</i> in terms of the <tt>splitCats</tt> that
375/// are listed.
376/// </p>
377///
378///
379/// <h4>Multiple connected builds</h4>
380/// <p>
381/// Sometimes you want to build multiple PDFs for independent consecutive fits
382/// that share some of their parameters. For example, we have two prototype PDFs
383/// <tt>pdfA(x;p,q)</tt> and <tt>pdfB(x;p,r)</tt> that have a common parameter <tt>p</tt>.
384/// We want to build a <tt>RooSimultaneous</tt> for both <tt>pdfA</tt> and <tt>B</tt>,
385/// which involves a split of parameter <tt>p</tt> and we would like to build the
386/// simultaneous pdfs <tt>simA</tt> and <tt>simB</tt> such that still share their (now split) parameters
387/// <tt>p_XXX</tt>. This is accomplished by letting a single instance of <tt>RooSimPdfBuilder</tt> handle
388/// the builds of both <tt>pdfA</tt> and <tt>pdfB</tt>, as illustrated in this example:
389/// </p>
390/// <pre>
391/// RooSimPdfBuilder builder(RooArgSet(pdfA,pdfB)) ;
392///
393/// RooArgSet* configA = builder.createProtoBuildConfig() ;
394/// (*configA)["physModels"] = "pdfA" ;
395/// (*configA)["splitCats"] = "C" ;
396/// (*configA)["pdf"] = "C : p" ;
397/// RooSimultaneous* simA = builder.buildPdf(*configA,&D) ;
398///
399/// RooArgSet* configB = builder.createProtoBuildConfig() ;
400/// (*configA)["physModels"] = "pdfB" ;
401/// (*configA)["splitCats"] = "C" ;
402/// (*configA)["pdf"] = "C : p" ;
403/// RooSimultaneous* simB = builder.buildPdf(*configB,&D) ;
404/// </pre>
405///
406/// <h2>Ownership of constructed PDFs</h2>
407/// <p>
408/// The <tt>RooSimPdfBuilder</tt> instance owns all the objects it creates, including the top-level
409/// <tt>RooSimultaneous</tt> returned by <tt>buildPdf()</tt>. Therefore the builder instance should
410/// exist as long as the constructed PDFs needs to exist.
411/// </p>
412///
413///
414///
415
416
417#include <cstring>
418#include "strtok.h"
419
420#ifndef _WIN32
421#include <strings.h>
422#endif
423
424#include "Riostream.h"
425#include "RooSimPdfBuilder.h"
426
427#include "RooRealVar.h"
428#include "RooFormulaVar.h"
429#include "RooAbsCategory.h"
430#include "RooCategory.h"
431#include "RooStringVar.h"
432#include "RooMappedCategory.h"
433#include "RooRealIntegral.h"
434#include "RooDataSet.h"
435#include "RooArgSet.h"
436#include "RooPlot.h"
437#include "RooAddPdf.h"
438#include "RooLinearVar.h"
439#include "RooTruthModel.h"
440#include "RooAddModel.h"
441#include "RooProdPdf.h"
442#include "RooCustomizer.h"
443#include "RooThresholdCategory.h"
444#include "RooMultiCategory.h"
445#include "RooSuperCategory.h"
446#include "RooSimultaneous.h"
447#include "RooTrace.h"
448#include "RooFitResult.h"
449#include "RooDataHist.h"
450#include "RooGenericPdf.h"
451#include "RooMsgService.h"
452
453#include "ROOT/StringUtils.hxx"
454
455using namespace std ;
456
458;
459
460
461
462////////////////////////////////////////////////////////////////////////////////
463
465 _protoPdfSet(protoPdfSet)
466{
470}
471
472
473
474
475////////////////////////////////////////////////////////////////////////////////
476/// Make RooArgSet of configuration objects
477
479{
480 RooArgSet* buildConfig = new RooArgSet ;
481 buildConfig->addOwned(* new RooStringVar("physModels","List and mapping of physics models to include in build","",4096)) ;
482 buildConfig->addOwned(* new RooStringVar("splitCats","List of categories used for splitting","",1024)) ;
483
484 for(auto * proto : _protoPdfSet) {
485 buildConfig->addOwned(* new RooStringVar(proto->GetName(),proto->GetName(),"",4096)) ;
486 }
487
488 return buildConfig ;
489}
490
491
492
493////////////////////////////////////////////////////////////////////////////////
494
496{
497 _splitNodeList.add(specSet) ;
498}
499
500
501
502////////////////////////////////////////////////////////////////////////////////
503/// Initialize needed components
504
505RooSimultaneous* RooSimPdfBuilder::buildPdf(const RooArgSet& buildConfig, const RooArgSet& dependents,
506 const RooArgSet* auxSplitCats, bool verbose)
507{
508 const char* spaceChars = " \t" ;
509
510 // Retrieve physics index category
511 int buflen = strlen(((RooStringVar*)buildConfig.find("physModels"))->getVal())+1 ;
512 std::vector<char> bufContainer(buflen);
513 char *buf = bufContainer.data() ;
514
515 strlcpy(buf,((RooStringVar*)buildConfig.find("physModels"))->getVal(),buflen) ;
516 RooAbsCategoryLValue* physCat(0) ;
517 if (strstr(buf," : ")) {
518 const char* physCatName = strtok(buf,spaceChars) ;
519 physCat = dynamic_cast<RooAbsCategoryLValue*>(dependents.find(physCatName)) ;
520 if (!physCat) {
521 coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR physics index category " << physCatName
522 << " not found in dataset variables" << endl ;
523 return 0 ;
524 }
525 coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: category indexing physics model: " << physCatName << endl ;
526 }
527
528 // Create list of physics models to be built
529 char *physName ;
530 RooArgSet physModelSet ;
531 if (physCat) {
532 // Absorb colon token
533 strtok(0,spaceChars) ;
534 physName = strtok(0,spaceChars) ;
535 } else {
536 physName = strtok(buf,spaceChars) ;
537 }
538
539 if (!physName) {
540 coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR: No models specified, nothing to do!" << endl ;
541 return 0 ;
542 }
543
544 bool first{true} ;
545 RooArgSet stateMap ;
546 while(physName) {
547
548 char *stateName(0) ;
549
550 // physName may be <state>=<pdfName> or just <pdfName> is state and pdf have identical names
551 if (strchr(physName,'=')) {
552 // Must have a physics category for mapping to make sense
553 if (!physCat) {
554 coutW(ObjectHandling) << "RooSimPdfBuilder::buildPdf: WARNING: without physCat specification "
555 << "<physCatState>=<pdfProtoName> association is meaningless" << endl ;
556 }
557 stateName = physName ;
558 physName = strchr(stateName,'=') ;
559 if (physName) {
560 *(physName++) = 0 ;
561 }
562 } else {
563 stateName = physName ;
564 }
565
566 RooAbsPdf* physModel = (RooAbsPdf*) (physName ? _protoPdfSet.find(physName) : 0 );
567 if (!physModel) {
568 coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR requested physics model "
569 << (physName?physName:"(null)") << " is not defined" << endl ;
570 return 0 ;
571 }
572
573 // Check if state mapping has already been defined
574 if (stateMap.find(stateName)) {
575 coutW(InputArguments) << "RooSimPdfBuilder::buildPdf: WARNING: multiple PDFs specified for state "
576 << stateName << ", only first will be used" << endl ;
577 continue ;
578 }
579
580 // Add pdf to list of models to be processed
581 physModelSet.add(*physModel,true) ; // silence duplicate insertion warnings
582
583 // Store state->pdf mapping
584 stateMap.addOwned(* new RooStringVar(stateName,stateName,physName)) ;
585
586 // Continue with next mapping
587 physName = strtok(0,spaceChars) ;
588 if (first) {
589 first = false ;
590 } else if (physCat==0) {
591 coutW(InputArguments) << "RooSimPdfBuilder::buildPdf: WARNING: without physCat specification, only the first model will be used" << endl ;
592 break ;
593 }
594 }
595 coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: list of physics models " << physModelSet << endl ;
596
597
598
599 // Create list of dataset categories to be used in splitting
600 TList splitStateList ;
601 RooArgSet splitCatSet ;
602
603 buflen = strlen(((RooStringVar*)buildConfig.find("splitCats"))->getVal())+1 ;
604 bufContainer.resize(buflen);
605 buf = bufContainer.data(); // memory might have been reallocated after resize, so reset the data pointer
606
607 strlcpy(buf,((RooStringVar*)buildConfig.find("splitCats"))->getVal(),buflen) ;
608
609 char *catName = strtok(buf,spaceChars) ;
610 char *stateList(0) ;
611 while(catName) {
612
613 // Chop off optional list of selected states
614 char* tokenPtr(0) ;
615 if (strchr(catName,'(')) {
616
617 catName = R__STRTOK_R(catName,"(",&tokenPtr) ;
618 stateList = R__STRTOK_R(0,")",&tokenPtr) ;
619
620 } else {
621 stateList = 0 ;
622 }
623
624 RooCategory* splitCat = catName ? dynamic_cast<RooCategory*>(dependents.find(catName)) : 0 ;
625 if (!splitCat) {
626 coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR requested split category " << (catName?catName:"(null)")
627 << " is not a RooCategory in the dataset" << endl ;
628 return 0 ;
629 }
630 splitCatSet.add(*splitCat) ;
631
632 // Process optional state list
633 if (stateList) {
634 coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: splitting of category " << catName
635 << " restricted to states (" << stateList << ")" << endl ;
636
637 // Create list named after this splitCat holding its selected states
638 TList* slist = new TList ;
639 slist->SetName(catName) ;
640 splitStateList.Add(slist) ;
641
642 char* stateLabel = R__STRTOK_R(stateList,",",&tokenPtr) ;
643
644 while(stateLabel) {
645 // Lookup state label and require it exists
646 const RooCatType* type = splitCat->lookupType(stateLabel) ;
647 if (!type) {
648 coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR splitCat " << splitCat->GetName()
649 << " doesn't have a state named " << stateLabel << endl ;
650 splitStateList.Delete() ;
651 return 0 ;
652 }
653 slist->Add((TObject*)type) ;
654
655 stateLabel = R__STRTOK_R(0,",",&tokenPtr) ;
656 }
657 }
658
659 catName = strtok(0,spaceChars) ;
660 }
661 if (physCat) splitCatSet.add(*physCat) ;
662 RooSuperCategory masterSplitCat("masterSplitCat","Master splitting category",splitCatSet) ;
663
664 coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: list of splitting categories " << splitCatSet << endl ;
665
666 // Clone auxiliary split cats and attach to splitCatSet
667 RooArgSet auxSplitSet ;
668 if (auxSplitCats) {
669 // Deep clone auxililary split cats
670 if (!std::unique_ptr<RooArgSet>{auxSplitCats->snapshot(true)}) {
671 coutE(InputArguments) << "RooSimPdfBuilder::buildPdf(" << GetName() << ") Couldn't deep-clone set auxiliary splitcats, abort." << endl ;
672 return 0 ;
673 }
674
675 for (const auto arg : *auxSplitCats) {
676 // Find counterpart in cloned set
677 RooAbsArg* aux = auxSplitCats->find(arg->GetName()) ;
678
679 // Check that there is no fundamental splitCat in the dataset with the bane of the auxiliary split
680 if (splitCatSet.find(aux->GetName())) {
681 coutW(InputArguments) << "RooSimPdfBuilder::buildPdf: WARNING: dataset contains a fundamental splitting category " << endl
682 << " with the same name as an auxiliary split function (" << aux->GetName() << "). " << endl
683 << " Auxiliary split function will be ignored" << endl ;
684 continue ;
685 }
686
687 // Check that all servers of this aux cat are contained in splitCatSet
688 std::unique_ptr<RooArgSet> parSet{aux->getParameters(splitCatSet)} ;
689 if (parSet->getSize()>0) {
690 coutW(InputArguments) << "RooSimPdfBuilder::buildPdf: WARNING: ignoring auxiliary category " << aux->GetName()
691 << " because it has servers that are not listed in splitCatSet: " << *parSet << endl ;
692 continue ;
693 }
694
695 // Redirect servers to splitCatSet
696 aux->recursiveRedirectServers(splitCatSet) ;
697
698 // Add top level nodes to auxSplitSet
699 auxSplitSet.add(*aux) ;
700 }
701
702 coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: list of auxiliary splitting categories " << auxSplitSet << endl ;
703 }
704
705 TList customizerList;
706
707 // Loop over requested physics models and build components
708 for (const auto arg : physModelSet) {
709 auto physModel = static_cast<const RooAbsPdf*>(arg);
710
711 coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: processing physics model " << physModel->GetName() << endl ;
712
713 RooCustomizer* physCustomizer = new RooCustomizer(*physModel,masterSplitCat,_splitNodeList) ;
714 customizerList.Add(physCustomizer) ;
715
716 // Parse the splitting rules for this physics model
717 RooStringVar* ruleStr = (RooStringVar*) buildConfig.find(physModel->GetName()) ;
718 if (ruleStr) {
719
720 enum Mode { SplitCat, Colon, ParamList } ;
721 Mode mode(SplitCat) ;
722
723 const char* splitCatName ;
724 RooAbsCategory* splitCat(0) ;
725
726 for (const auto& token : ROOT::Split(ruleStr->getVal(), spaceChars)) {
727
728 switch (mode) {
729 case SplitCat:
730 {
731 splitCatName = token.data();
732
733 if (token.find(',') != std::string::npos) {
734 // Composite splitting category
735
736 // Check if already instantiated
737 splitCat = (RooAbsCategory*) _compSplitCatSet.find(splitCatName) ;
738 const std::string origCompCatName{splitCatName} ;
739 if (!splitCat) {
740 // Build now
741
742 RooArgSet compCatSet ;
743 for (const auto& catName2 : ROOT::Split(token, ",")) {
744 RooAbsArg* cat = splitCatSet.find(catName2.data()) ;
745
746 // If not, check if it is an auxiliary splitcat
747 if (!cat) {
748 cat = (RooAbsCategory*) auxSplitSet.find(catName2.data()) ;
749 }
750
751 if (!cat) {
752 coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR " << catName2
753 << " not found in the primary or auxilary splitcat list" << endl ;
754 customizerList.Delete() ;
755
756 splitStateList.Delete() ;
757 return 0 ;
758 }
759 compCatSet.add(*cat) ;
760 }
761
762
763 // check that any auxSplitCats in compCatSet do not depend on any other
764 // fundamental or auxiliary splitting categories in the composite set.
765 for (const auto theArg : compCatSet) {
766 RooArgSet tmp(compCatSet) ;
767 tmp.remove(*theArg) ;
768 if (theArg->dependsOnValue(tmp)) {
770 << "RooSimPdfBuilder::buildPDF: ERROR: Ill defined split: auxiliary splitting category "
771 << theArg->GetName() << " used in composite split " << compCatSet
772 << " depends on one or more of the other splitting categories in the composite split"
773 << std::endl ;
774
775 // Cleanup and axit
776 customizerList.Delete() ;
777 splitStateList.Delete() ;
778 return 0 ;
779 }
780 }
781
782 splitCat = new RooMultiCategory(origCompCatName.c_str(), origCompCatName.c_str(), compCatSet) ;
783 _compSplitCatSet.addOwned(*splitCat) ;
784 //cout << "composite splitcat: " << splitCat->GetName() ;
785 }
786 } else {
787 // Simple splitting category
788
789 // First see if it is a simple splitting category
790 splitCat = (RooAbsCategory*) splitCatSet.find(splitCatName) ;
791
792 // If not, check if it is an auxiliary splitcat
793 if (!splitCat) {
794 splitCat = (RooAbsCategory*) auxSplitSet.find(splitCatName) ;
795 }
796
797 if (!splitCat) {
798 coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR splitting category "
799 << splitCatName << " not found in the primary or auxiliary splitcat list" << endl ;
800 customizerList.Delete() ;
801 splitStateList.Delete() ;
802 return 0 ;
803 }
804 }
805
806 mode = Colon ;
807 break ;
808 }
809 case Colon:
810 {
811 if (token != ":") {
812 coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR in parsing, expected ':' after "
813 << splitCat << ", found " << token << endl ;
814 customizerList.Delete() ;
815 splitStateList.Delete() ;
816 return 0 ;
817 }
818 mode = ParamList ;
819 break ;
820 }
821 case ParamList:
822 {
823 // Verify the validity of the parameter list and build the corresponding argset
824 RooArgSet splitParamList ;
825 std::unique_ptr<RooArgSet> paramList{physModel->getParameters(dependents)} ;
826
827 // wve -- add nodes to parameter list
828 std::unique_ptr<RooArgSet> compList{physModel->getComponents()} ;
829 paramList->add(*compList) ;
830
831 const bool lastCharIsComma = (token[token.size()-1]==',') ;
832
833 for (const auto& paramName : ROOT::Split(token, ",", /*skipEmpty= */ true)) {
834 // Check for fractional split option 'param_name[remainder_state]'
835 std::string remainderState;
836 {
837 const auto pos = paramName.find('[');
838 if (pos != std::string::npos) {
839 const auto posEnd = paramName.find(']');
840 remainderState = paramName.substr(pos+1, posEnd - pos - 1);
841 }
842 }
843
844
845 // If fractional split is specified, check that remainder state is a valid state of this split cat
846 if (!remainderState.empty()) {
847 if (!splitCat->hasLabel(remainderState)) {
848 coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR fraction split of parameter "
849 << paramName << " has invalid remainder state name: " << remainderState << endl ;
850 customizerList.Delete() ;
851 splitStateList.Delete() ;
852 return 0 ;
853 }
854 }
855
856 RooAbsArg* param = paramList->find(paramName.data()) ;
857 if (!param) {
858 coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR " << paramName
859 << " is not a parameter of physics model " << physModel->GetName() << endl ;
860 customizerList.Delete() ;
861 splitStateList.Delete() ;
862 return 0 ;
863 }
864 splitParamList.add(*param) ;
865
866 // Build split leaf of fraction splits here
867 if (!remainderState.empty()) {
868
869 // Check if we are splitting a real-valued parameter
870 if (!dynamic_cast<RooAbsReal*>(param)) {
871 coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR fraction split requested of non-real valued parameter "
872 << param->GetName() << endl ;
873 customizerList.Delete() ;
874 splitStateList.Delete() ;
875 return 0 ;
876 }
877
878 // Check if we are doing a restricted build
879 TList* remStateSplitList = static_cast<TList*>(splitStateList.FindObject(splitCat->GetName())) ;
880
881 // If so, check if remainder state is actually being built.
882 if (remStateSplitList && !remStateSplitList->FindObject(remainderState.data())) {
883 coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR " << paramName
884 << " remainder state " << remainderState << " in parameter split "
885 << param->GetName() << " is not actually being built" << endl ;
886 customizerList.Delete() ;
887 splitStateList.Delete() ;
888 return 0 ;
889 }
890
891 std::unique_ptr<TIterator> iter{splitCat->typeIterator()} ;
893 RooArgList fracLeafList ;
894 std::string formExpr{"1"} ;
895 int i(0) ;
896
897 while((type=(RooCatType*)iter->Next())) {
898
899 // Skip remainder state
900 if (std::string{type->GetName()} != remainderState) continue ;
901
902 // If restricted build is requested, skip states of splitcat that are not built
903 if (remStateSplitList && !remStateSplitList->FindObject(type->GetName())) {
904 continue ;
905 }
906
907 // Construct name of split leaf
908 const auto splitLeafName = std::string{param->GetName()} + "_" + type->GetName();
909
910 // Check if split leaf already exists
911 RooAbsArg* splitLeaf = _splitNodeList.find(splitLeafName.c_str()) ;
912 if (!splitLeaf) {
913 // If not create it now
914 splitLeaf = (RooAbsArg*) param->clone(splitLeafName.c_str()) ;
915 _splitNodeList.add(*splitLeaf) ;
916 _splitNodeListOwned.addOwned(*splitLeaf) ;
917 }
918 fracLeafList.add(*splitLeaf) ;
919 formExpr += Form("-@%d",i++) ;
920 }
921
922 // Construct RooFormulaVar expresssing remainder of fraction
923 const auto remLeafName = std::string{param->GetName()} + "_" + remainderState;
924
925 // Check if no specialization was already specified for remainder state
926 if (!_splitNodeList.find(remLeafName.c_str())) {
927 RooAbsArg* remLeaf = new RooFormulaVar(remLeafName.c_str(), formExpr.c_str(), fracLeafList) ;
928 _splitNodeList.add(*remLeaf) ;
929 _splitNodeListOwned.addOwned(*remLeaf) ;
930 coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: creating remainder fraction formula for " << remainderState
931 << " specialization of split parameter " << param->GetName() << " " << formExpr << endl ;
932 }
933 }
934 }
935
936 // Add the rule to the appropriate customizer ;
937 physCustomizer->splitArgs(splitParamList,*splitCat) ;
938
939 if (!lastCharIsComma) mode = SplitCat ;
940 break ;
941 }
942 }
943 }
944 if (mode!=SplitCat) {
945 coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR in parsing, expected "
946 << (mode==Colon?":":"parameter list") << " in " << ruleStr->getVal() << endl ;
947 }
948
949 //RooArgSet* paramSet = physModel->getParameters(dependents) ;
950 } else {
951 coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: no splitting rules for " << physModel->GetName() << endl ;
952 }
953 }
954
955 coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: configured customizers for all physics models" << endl ;
957 customizerList.Print() ;
958 }
959
960 // Create fit category from physCat and splitCatList ;
961 RooArgSet fitCatList ;
962 if (physCat) fitCatList.add(*physCat) ;
963 fitCatList.add(splitCatSet) ;
964 RooSuperCategory *fitCat = new RooSuperCategory("fitCat","fitCat",fitCatList) ;
965
966 // Create master PDF
967 RooSimultaneous* simPdf = new RooSimultaneous("simPdf","simPdf",*fitCat) ;
968
969 // Add component PDFs to master PDF
970 std::unique_ptr<TIterator> fcIter{fitCat->typeIterator()} ;
971
972 RooCatType* fcState ;
973 while((fcState=(RooCatType*)fcIter->Next())) {
974 // Select fitCat state
975 fitCat->setLabel(fcState->GetName()) ;
976
977 // Check if this fitCat state is selected
978 bool select{true} ;
979 for (const auto arg : fitCatList) {
980 auto splitCat = static_cast<const RooAbsCategory*>(arg);
981
982 // Find selected state list
983 TList* slist = (TList*) splitStateList.FindObject(splitCat->GetName()) ;
984 if (!slist) continue ;
985 RooCatType* type = (RooCatType*) slist->FindObject(splitCat->getCurrentLabel()) ;
986 if (!type) {
987 select = false ;
988 }
989 }
990 if (!select) continue ;
991
992
993 // Select appropriate PDF for this physCat state
994 RooCustomizer* physCustomizer ;
995 if (physCat) {
996 RooStringVar* physNameVar = (RooStringVar*) stateMap.find(physCat->getCurrentLabel()) ;
997 if (!physNameVar) continue ;
998 physCustomizer = static_cast<RooCustomizer*>(customizerList.FindObject(physNameVar->getVal())) ;
999 } else {
1000 physCustomizer = static_cast<RooCustomizer*>(customizerList.First()) ;
1001 }
1002
1003 coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: Customizing physics model " << physCustomizer->GetName()
1004 << " for mode " << fcState->GetName() << endl ;
1005
1006 // Customizer PDF for current state and add to master simPdf
1007 RooAbsPdf* fcPdf = (RooAbsPdf*) physCustomizer->build(masterSplitCat.getCurrentLabel(),verbose) ;
1008 simPdf->addPdf(*fcPdf,fcState->GetName()) ;
1009 }
1010
1011 // Move customizers (owning the cloned branch node components) to the attic
1012 _retiredCustomizerList.AddAll(&customizerList) ;
1013
1014 splitStateList.Delete() ;
1015
1016 _simPdfList.push_back(simPdf) ;
1017 _fitCatList.push_back(fitCat) ;
1018 return simPdf ;
1019}
1020
1021
1022
1023
1024
1025////////////////////////////////////////////////////////////////////////////////
1026
1028{
1030
1031 for(auto * ptr : _simPdfList) {
1032 delete ptr;
1033 }
1034 for(auto * ptr : _fitCatList) {
1035 delete ptr;
1036 }
1037}
1038
1039
#define coutI(a)
Definition: RooMsgService.h:30
#define coutW(a)
Definition: RooMsgService.h:32
#define coutE(a)
Definition: RooMsgService.h:33
#define oodologI(o, a)
Definition: RooMsgService.h:73
#define ClassImp(name)
Definition: Rtypes.h:375
Option_t Option_t TPoint TPoint const char mode
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 * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition: TString.cxx:2447
const char * proto
Definition: civetweb.c:17493
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:78
RooArgSet * getParameters(const RooAbsData *data, bool stripDisconnected=true) const
Create a list of leaf nodes in the arg tree starting with ourself as top node that don't match any of...
Definition: RooAbsArg.cxx:569
Bool_t recursiveRedirectServers(const RooAbsCollection &newServerList, Bool_t mustReplaceAll=kFALSE, Bool_t nameChange=kFALSE, Bool_t recurseInNewSet=kTRUE)
Recursively replace all servers with the new servers in newSet.
Definition: RooAbsArg.cxx:1197
virtual TObject * clone(const char *newname=0) const =0
RooAbsCategoryLValue is the common abstract base class for objects that represent a discrete value th...
A space to attach TBranches.
bool hasLabel(const std::string &label) const
Check if a state with name label exists.
virtual const char * getCurrentLabel() const
Return label string of current state.
TIterator * typeIterator() const
const RooCatType * lookupType(value_type index, Bool_t printError=kFALSE) const
Find our type corresponding to the specified index, or return nullptr for no match.
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
virtual Bool_t addOwned(RooAbsArg &var, Bool_t silent=kFALSE)
Add an argument and transfer the ownership to the collection.
void useHashMapForFind(bool flag) const
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE)
Remove the specified argument from our list.
RooAbsArg * find(const char *name) const
Find object with given name in list.
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:63
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:22
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:57
RooArgSet * snapshot(bool deepCopy=true) const
Use RooAbsCollection::snapshot(), but return as RooArgSet.
Definition: RooArgSet.h:180
RooCatType is an auxilary class for RooAbsCategory and defines a a single category state.
const Text_t * GetName() const override
Returns name of object.
RooCategory is an object to represent discrete states.
Definition: RooCategory.h:28
RooCustomizer is a factory class to produce clones of a prototype composite PDF object with the same ...
Definition: RooCustomizer.h:35
void splitArgs(const RooArgSet &argSet, const RooAbsCategory &splitCat)
Split all arguments in 'set' into individualized clones for each defined state of 'splitCat'.
RooAbsArg * build(const char *masterCatState, Bool_t verbose=kFALSE)
Build a clone of the prototype executing all registered 'replace' rules and 'split' rules for the mas...
A RooFormulaVar is a generic implementation of a real-valued object, which takes a RooArgList of serv...
Definition: RooFormulaVar.h:30
RooMultiCategory connects several RooAbsCategory objects into a single category.
This tool has now been superseded by RooSimWSTool
RooArgSet _compSplitCatSet
List of owned composite splitting categories.
RooArgSet _protoPdfSet
Set of prototype PDFS.
std::list< RooSuperCategory * > _fitCatList
The super-categories that we built.
RooArgSet _splitNodeList
List of owned split nodes.
std::list< RooSimultaneous * > _simPdfList
The simpdfs that we built.
RooArgSet * createProtoBuildConfig()
Make RooArgSet of configuration objects.
void addSpecializations(const RooArgSet &specSet)
~RooSimPdfBuilder() override
RooArgSet _splitNodeListOwned
List of all split nodes.
TList _retiredCustomizerList
Retired customizer from previous builds (own their PDF branch nodes)
RooSimultaneous * buildPdf(const RooArgSet &buildConfig, const RooArgSet &dependents, const RooArgSet *auxSplitCats=0, bool verbose=false)
Initialize needed components.
RooSimPdfBuilder(const RooArgSet &pdfProtoList)
RooSimultaneous facilitates simultaneous fitting of multiple PDFs to subsets of a given dataset.
Bool_t addPdf(const RooAbsPdf &pdf, const char *catLabel)
Associate given PDF with index category state label 'catLabel'.
RooStringVar is a RooAbsArg implementing string values.
Definition: RooStringVar.h:23
const char * getVal() const
Definition: RooStringVar.h:34
The RooSuperCategory can join several RooAbsCategoryLValue objects into a single category.
Bool_t setLabel(const char *label, Bool_t printError=kTRUE) override
Set the value of the super category by specifying the state name.
void SetName(const char *name)
Definition: TCollection.h:206
virtual void AddAll(const TCollection *col)
Add all objects from collection col to this collection.
void Print(Option_t *option="") const override
Default print for collections, calls Print(option, 1).
A doubly linked list.
Definition: TList.h:44
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:87
TObject * First() const override
Return the first object in the list. Returns 0 when list is empty.
Definition: TList.cxx:659
void Delete(Option_t *option="") override
Remove all objects from the list AND delete all heap based objects.
Definition: TList.cxx:470
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
Mother of all ROOT objects.
Definition: TObject.h:37
virtual const char * GetName() const
Returns name of object.
Definition: TObject.cxx:359
std::vector< std::string > Split(std::string_view str, std::string_view delims, bool skipEmpty=false)
Splits a string at each character in delims.
Definition: StringUtils.cxx:23
@ InputArguments
Definition: RooGlobalFunc.h:64
@ ObjectHandling
Definition: RooGlobalFunc.h:64
Definition: first.py:1