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* config = builder.createProtoBuildConfig() ;
88 /// (*config)["physModels"] = "pdf" ; /// Name of the PDF we are going to work with
89 /// (*config)["splitCats"] = "C" ; /// Category used to differentiate sub-datasets
90 /// (*config)["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(*config,&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_t verbose=kFALSE) {
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 "RooFit.h"
418 
419 #include <cstring>
420 #include "strtok.h"
421 
422 #ifndef _WIN32
423 #include <strings.h>
424 #endif
425 
426 #include "Riostream.h"
427 #include "RooSimPdfBuilder.h"
428 
429 #include "RooRealVar.h"
430 #include "RooFormulaVar.h"
431 #include "RooAbsCategory.h"
432 #include "RooCategory.h"
433 #include "RooStringVar.h"
434 #include "RooMappedCategory.h"
435 #include "RooRealIntegral.h"
436 #include "RooDataSet.h"
437 #include "RooArgSet.h"
438 #include "RooPlot.h"
439 #include "RooAddPdf.h"
440 #include "RooLinearVar.h"
441 #include "RooTruthModel.h"
442 #include "RooAddModel.h"
443 #include "RooProdPdf.h"
444 #include "RooCustomizer.h"
445 #include "RooThresholdCategory.h"
446 #include "RooMultiCategory.h"
447 #include "RooSuperCategory.h"
448 #include "RooSimultaneous.h"
449 #include "RooTrace.h"
450 #include "RooFitResult.h"
451 #include "RooDataHist.h"
452 #include "RooGenericPdf.h"
453 #include "RooMsgService.h"
454 #include "RooHelpers.h"
455 
456 using namespace std ;
457 
459 ;
460 
461 
462 
463 ////////////////////////////////////////////////////////////////////////////////
464 
466  _protoPdfSet(protoPdfSet)
467 {
471 }
472 
473 
474 
475 
476 ////////////////////////////////////////////////////////////////////////////////
477 /// Make RooArgSet of configuration objects
478 
480 {
481  RooArgSet* buildConfig = new RooArgSet ;
482  buildConfig->addOwned(* new RooStringVar("physModels","List and mapping of physics models to include in build","",4096)) ;
483  buildConfig->addOwned(* new RooStringVar("splitCats","List of categories used for splitting","",1024)) ;
484 
486  RooAbsPdf* proto ;
487  while ((proto=(RooAbsPdf*)iter->Next())) {
488  buildConfig->addOwned(* new RooStringVar(proto->GetName(),proto->GetName(),"",4096)) ;
489  }
490  delete iter ;
491 
492  return buildConfig ;
493 }
494 
495 
496 
497 ////////////////////////////////////////////////////////////////////////////////
498 
500 {
501  _splitNodeList.add(specSet) ;
502 }
503 
504 
505 
506 ////////////////////////////////////////////////////////////////////////////////
507 /// Initialize needed components
508 
509 RooSimultaneous* RooSimPdfBuilder::buildPdf(const RooArgSet& buildConfig, const RooArgSet& dependents,
510  const RooArgSet* auxSplitCats, Bool_t verbose)
511 {
512  const char* spaceChars = " \t" ;
513 
514  // Retrieve physics index category
515  Int_t buflen = strlen(((RooStringVar*)buildConfig.find("physModels"))->getVal())+1 ;
516  char *buf = new char[buflen] ;
517 
518  strlcpy(buf,((RooStringVar*)buildConfig.find("physModels"))->getVal(),buflen) ;
519  RooAbsCategoryLValue* physCat(0) ;
520  if (strstr(buf," : ")) {
521  const char* physCatName = strtok(buf,spaceChars) ;
522  physCat = dynamic_cast<RooAbsCategoryLValue*>(dependents.find(physCatName)) ;
523  if (!physCat) {
524  coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR physics index category " << physCatName
525  << " not found in dataset variables" << endl ;
526  delete[] buf ;
527  return 0 ;
528  }
529  coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: category indexing physics model: " << physCatName << endl ;
530  }
531 
532  // Create list of physics models to be built
533  char *physName ;
534  RooArgSet physModelSet ;
535  if (physCat) {
536  // Absorb colon token
537  strtok(0,spaceChars) ;
538  physName = strtok(0,spaceChars) ;
539  } else {
540  physName = strtok(buf,spaceChars) ;
541  }
542 
543  if (!physName) {
544  coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR: No models specified, nothing to do!" << endl ;
545  delete[] buf ;
546  return 0 ;
547  }
548 
549  Bool_t first(kTRUE) ;
550  RooArgSet stateMap ;
551  while(physName) {
552 
553  char *stateName(0) ;
554 
555  // physName may be <state>=<pdfName> or just <pdfName> is state and pdf have identical names
556  if (strchr(physName,'=')) {
557  // Must have a physics category for mapping to make sense
558  if (!physCat) {
559  coutW(ObjectHandling) << "RooSimPdfBuilder::buildPdf: WARNING: without physCat specification "
560  << "<physCatState>=<pdfProtoName> association is meaningless" << endl ;
561  }
562  stateName = physName ;
563  physName = strchr(stateName,'=') ;
564  if (physName) {
565  *(physName++) = 0 ;
566  }
567  } else {
568  stateName = physName ;
569  }
570 
571  RooAbsPdf* physModel = (RooAbsPdf*) (physName ? _protoPdfSet.find(physName) : 0 );
572  if (!physModel) {
573  coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR requested physics model "
574  << (physName?physName:"(null)") << " is not defined" << endl ;
575  delete[] buf ;
576  return 0 ;
577  }
578 
579  // Check if state mapping has already been defined
580  if (stateMap.find(stateName)) {
581  coutW(InputArguments) << "RooSimPdfBuilder::buildPdf: WARNING: multiple PDFs specified for state "
582  << stateName << ", only first will be used" << endl ;
583  continue ;
584  }
585 
586  // Add pdf to list of models to be processed
587  physModelSet.add(*physModel,kTRUE) ; // silence duplicate insertion warnings
588 
589  // Store state->pdf mapping
590  stateMap.addOwned(* new RooStringVar(stateName,stateName,physName)) ;
591 
592  // Continue with next mapping
593  physName = strtok(0,spaceChars) ;
594  if (first) {
595  first = kFALSE ;
596  } else if (physCat==0) {
597  coutW(InputArguments) << "RooSimPdfBuilder::buildPdf: WARNING: without physCat specification, only the first model will be used" << endl ;
598  break ;
599  }
600  }
601  coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: list of physics models " << physModelSet << endl ;
602 
603 
604 
605  // Create list of dataset categories to be used in splitting
606  TList splitStateList ;
607  RooArgSet splitCatSet ;
608 
609  delete[] buf ;
610  buflen = strlen(((RooStringVar*)buildConfig.find("splitCats"))->getVal())+1 ;
611  buf = new char[buflen] ;
612  strlcpy(buf,((RooStringVar*)buildConfig.find("splitCats"))->getVal(),buflen) ;
613 
614  char *catName = strtok(buf,spaceChars) ;
615  char *stateList(0) ;
616  while(catName) {
617 
618  // Chop off optional list of selected states
619  char* tokenPtr(0) ;
620  if (strchr(catName,'(')) {
621 
622  catName = R__STRTOK_R(catName,"(",&tokenPtr) ;
623  stateList = R__STRTOK_R(0,")",&tokenPtr) ;
624 
625  } else {
626  stateList = 0 ;
627  }
628 
629  RooCategory* splitCat = catName ? dynamic_cast<RooCategory*>(dependents.find(catName)) : 0 ;
630  if (!splitCat) {
631  coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR requested split category " << (catName?catName:"(null)")
632  << " is not a RooCategory in the dataset" << endl ;
633  delete[] buf ;
634  return 0 ;
635  }
636  splitCatSet.add(*splitCat) ;
637 
638  // Process optional state list
639  if (stateList) {
640  coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: splitting of category " << catName
641  << " restricted to states (" << stateList << ")" << endl ;
642 
643  // Create list named after this splitCat holding its selected states
644  TList* slist = new TList ;
645  slist->SetName(catName) ;
646  splitStateList.Add(slist) ;
647 
648  char* stateLabel = R__STRTOK_R(stateList,",",&tokenPtr) ;
649 
650  while(stateLabel) {
651  // Lookup state label and require it exists
652  const RooCatType* type = splitCat->lookupType(stateLabel) ;
653  if (!type) {
654  coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR splitCat " << splitCat->GetName()
655  << " doesn't have a state named " << stateLabel << endl ;
656  splitStateList.Delete() ;
657  delete[] buf ;
658  return 0 ;
659  }
660  slist->Add((TObject*)type) ;
661 
662  stateLabel = R__STRTOK_R(0,",",&tokenPtr) ;
663  }
664  }
665 
666  catName = strtok(0,spaceChars) ;
667  }
668  if (physCat) splitCatSet.add(*physCat) ;
669  RooSuperCategory masterSplitCat("masterSplitCat","Master splitting category",splitCatSet) ;
670 
671  coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: list of splitting categories " << splitCatSet << endl ;
672 
673  // Clone auxiliary split cats and attach to splitCatSet
674  RooArgSet auxSplitSet ;
675  RooArgSet* auxSplitCloneSet(0) ;
676  if (auxSplitCats) {
677  // Deep clone auxililary split cats
678  auxSplitCloneSet = (RooArgSet*) auxSplitCats->snapshot(kTRUE) ;
679  if (!auxSplitCloneSet) {
680  coutE(InputArguments) << "RooSimPdfBuilder::buildPdf(" << GetName() << ") Couldn't deep-clone set auxiliary splitcats, abort." << endl ;
681  delete[] buf ;
682  return 0 ;
683  }
684 
685  for (const auto arg : *auxSplitCats) {
686  // Find counterpart in cloned set
687  RooAbsArg* aux = auxSplitCats->find(arg->GetName()) ;
688 
689  // Check that there is no fundamental splitCat in the dataset with the bane of the auxiliary split
690  if (splitCatSet.find(aux->GetName())) {
691  coutW(InputArguments) << "RooSimPdfBuilder::buildPdf: WARNING: dataset contains a fundamental splitting category " << endl
692  << " with the same name as an auxiliary split function (" << aux->GetName() << "). " << endl
693  << " Auxiliary split function will be ignored" << endl ;
694  continue ;
695  }
696 
697  // Check that all servers of this aux cat are contained in splitCatSet
698  RooArgSet* parSet = aux->getParameters(splitCatSet) ;
699  if (parSet->getSize()>0) {
700  coutW(InputArguments) << "RooSimPdfBuilder::buildPdf: WARNING: ignoring auxiliary category " << aux->GetName()
701  << " because it has servers that are not listed in splitCatSet: " << *parSet << endl ;
702  delete parSet ;
703  continue ;
704  }
705 
706  // Redirect servers to splitCatSet
707  aux->recursiveRedirectServers(splitCatSet) ;
708 
709  // Add top level nodes to auxSplitSet
710  auxSplitSet.add(*aux) ;
711  }
712 
713  coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: list of auxiliary splitting categories " << auxSplitSet << endl ;
714  }
715 
716  delete[] buf;
717 
718  TList* customizerList = new TList ;
719 
720  // Loop over requested physics models and build components
721  for (const auto arg : physModelSet) {
722  auto physModel = static_cast<const RooAbsPdf*>(arg);
723 
724  coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: processing physics model " << physModel->GetName() << endl ;
725 
726  RooCustomizer* physCustomizer = new RooCustomizer(*physModel,masterSplitCat,_splitNodeList) ;
727  customizerList->Add(physCustomizer) ;
728 
729  // Parse the splitting rules for this physics model
730  RooStringVar* ruleStr = (RooStringVar*) buildConfig.find(physModel->GetName()) ;
731  if (ruleStr) {
732 
733  enum Mode { SplitCat, Colon, ParamList } ;
734  Mode mode(SplitCat) ;
735 
736  const char* splitCatName ;
737  RooAbsCategory* splitCat(0) ;
738 
739  for (const auto& token : RooHelpers::tokenise(ruleStr->getVal(), spaceChars)) {
740 
741  switch (mode) {
742  case SplitCat:
743  {
744  splitCatName = token.data();
745 
746  if (token.find(',') != std::string::npos) {
747  // Composite splitting category
748 
749  // Check if already instantiated
750  splitCat = (RooAbsCategory*) _compSplitCatSet.find(splitCatName) ;
751  TString origCompCatName(splitCatName) ;
752  if (!splitCat) {
753  // Build now
754 
755  RooArgSet compCatSet ;
756  for (const auto& catName2 : RooHelpers::tokenise(token, ",")) {
757  RooAbsArg* cat = splitCatSet.find(catName2.data()) ;
758 
759  // If not, check if it is an auxiliary splitcat
760  if (!cat) {
761  cat = (RooAbsCategory*) auxSplitSet.find(catName2.data()) ;
762  }
763 
764  if (!cat) {
765  coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR " << catName2
766  << " not found in the primary or auxilary splitcat list" << endl ;
767  customizerList->Delete() ;
768  delete customizerList ;
769 
770  splitStateList.Delete() ;
771  return 0 ;
772  }
773  compCatSet.add(*cat) ;
774  }
775 
776 
777  // check that any auxSplitCats in compCatSet do not depend on any other
778  // fundamental or auxiliary splitting categories in the composite set.
779  for (const auto theArg : compCatSet) {
780  RooArgSet tmp(compCatSet) ;
781  tmp.remove(*theArg) ;
782  if (theArg->dependsOnValue(tmp)) {
783  coutE(InputArguments) << "RooSimPdfBuilder::buildPDF: ERROR: Ill defined split: auxiliary splitting category " << theArg->GetName()
784  << " used in composite split " << compCatSet << " depends on one or more of the other splitting categories in the composite split" << endl ;
785 
786  // Cleanup and axit
787  customizerList->Delete() ;
788  delete customizerList ;
789  splitStateList.Delete() ;
790  return 0 ;
791  }
792  }
793 
794  splitCat = new RooMultiCategory(origCompCatName,origCompCatName,compCatSet) ;
795  _compSplitCatSet.addOwned(*splitCat) ;
796  //cout << "composite splitcat: " << splitCat->GetName() ;
797  }
798  } else {
799  // Simple splitting category
800 
801  // First see if it is a simple splitting category
802  splitCat = (RooAbsCategory*) splitCatSet.find(splitCatName) ;
803 
804  // If not, check if it is an auxiliary splitcat
805  if (!splitCat) {
806  splitCat = (RooAbsCategory*) auxSplitSet.find(splitCatName) ;
807  }
808 
809  if (!splitCat) {
810  coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR splitting category "
811  << splitCatName << " not found in the primary or auxiliary splitcat list" << endl ;
812  customizerList->Delete() ;
813  delete customizerList ;
814  splitStateList.Delete() ;
815  return 0 ;
816  }
817  }
818 
819  mode = Colon ;
820  break ;
821  }
822  case Colon:
823  {
824  if (token != ":") {
825  coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR in parsing, expected ':' after "
826  << splitCat << ", found " << token << endl ;
827  customizerList->Delete() ;
828  delete customizerList ;
829  splitStateList.Delete() ;
830  return 0 ;
831  }
832  mode = ParamList ;
833  break ;
834  }
835  case ParamList:
836  {
837  // Verify the validity of the parameter list and build the corresponding argset
838  RooArgSet splitParamList ;
839  RooArgSet* paramList = physModel->getParameters(dependents) ;
840 
841  // wve -- add nodes to parameter list
842  RooArgSet* compList = physModel->getComponents() ;
843  paramList->add(*compList) ;
844  delete compList ;
845 
846  const bool lastCharIsComma = (token[token.size()-1]==',') ;
847 
848  for (const auto& paramName : RooHelpers::tokenise(token, ",")) {
849  // Check for fractional split option 'param_name[remainder_state]'
850  std::string remainderState;
851  if (const auto pos = paramName.find('[') != std::string::npos) {
852  const auto posEnd = paramName.find(']');
853  remainderState = paramName.substr(pos+1, posEnd - pos - 1);
854  }
855 
856 
857  // If fractional split is specified, check that remainder state is a valid state of this split cat
858  if (!remainderState.empty()) {
859  if (!splitCat->hasLabel(remainderState)) {
860  coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR fraction split of parameter "
861  << paramName << " has invalid remainder state name: " << remainderState << endl ;
862  delete paramList ;
863  customizerList->Delete() ;
864  delete customizerList ;
865  splitStateList.Delete() ;
866  return 0 ;
867  }
868  }
869 
870  RooAbsArg* param = paramList->find(paramName.data()) ;
871  if (!param) {
872  coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR " << paramName
873  << " is not a parameter of physics model " << physModel->GetName() << endl ;
874  delete paramList ;
875  customizerList->Delete() ;
876  delete customizerList ;
877  splitStateList.Delete() ;
878  return 0 ;
879  }
880  splitParamList.add(*param) ;
881 
882  // Build split leaf of fraction splits here
883  if (!remainderState.empty()) {
884 
885  // Check if we are splitting a real-valued parameter
886  if (!dynamic_cast<RooAbsReal*>(param)) {
887  coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR fraction split requested of non-real valued parameter "
888  << param->GetName() << endl ;
889  delete paramList ;
890  customizerList->Delete() ;
891  delete customizerList ;
892  splitStateList.Delete() ;
893  return 0 ;
894  }
895 
896  // Check if we are doing a restricted build
897  TList* remStateSplitList = static_cast<TList*>(splitStateList.FindObject(splitCat->GetName())) ;
898 
899  // If so, check if remainder state is actually being built.
900  if (remStateSplitList && !remStateSplitList->FindObject(remainderState.data())) {
901  coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR " << paramName
902  << " remainder state " << remainderState << " in parameter split "
903  << param->GetName() << " is not actually being built" << endl ;
904  delete paramList ;
905  customizerList->Delete() ;
906  delete customizerList ;
907  splitStateList.Delete() ;
908  return 0 ;
909  }
910 
911  TIterator* iter = splitCat->typeIterator() ;
912  RooCatType* type ;
913  RooArgList fracLeafList ;
914  TString formExpr("1") ;
915  Int_t i(0) ;
916 
917  while((type=(RooCatType*)iter->Next())) {
918 
919  // Skip remainder state
920  if (!TString(type->GetName()).CompareTo(remainderState)) continue ;
921 
922  // If restricted build is requested, skip states of splitcat that are not built
923  if (remStateSplitList && !remStateSplitList->FindObject(type->GetName())) {
924  continue ;
925  }
926 
927  // Construct name of split leaf
928  TString splitLeafName(param->GetName()) ;
929  splitLeafName.Append("_") ;
930  splitLeafName.Append(type->GetName()) ;
931 
932  // Check if split leaf already exists
933  RooAbsArg* splitLeaf = _splitNodeList.find(splitLeafName) ;
934  if (!splitLeaf) {
935  // If not create it now
936  splitLeaf = (RooAbsArg*) param->clone(splitLeafName) ;
937  _splitNodeList.add(*splitLeaf) ;
938  _splitNodeListOwned.addOwned(*splitLeaf) ;
939  }
940  fracLeafList.add(*splitLeaf) ;
941  formExpr.Append(Form("-@%d",i++)) ;
942  }
943  delete iter ;
944 
945  // Construct RooFormulaVar expresssing remainder of fraction
946  TString remLeafName(param->GetName()) ;
947  remLeafName.Append("_") ;
948  remLeafName.Append(remainderState) ;
949 
950  // Check if no specialization was already specified for remainder state
951  if (!_splitNodeList.find(remLeafName)) {
952  RooAbsArg* remLeaf = new RooFormulaVar(remLeafName,formExpr,fracLeafList) ;
953  _splitNodeList.add(*remLeaf) ;
954  _splitNodeListOwned.addOwned(*remLeaf) ;
955  coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: creating remainder fraction formula for " << remainderState
956  << " specialization of split parameter " << param->GetName() << " " << formExpr << endl ;
957  }
958  }
959  }
960 
961  // Add the rule to the appropriate customizer ;
962  physCustomizer->splitArgs(splitParamList,*splitCat) ;
963 
964  delete paramList ;
965 
966  if (!lastCharIsComma) mode = SplitCat ;
967  break ;
968  }
969  }
970  }
971  if (mode!=SplitCat) {
972  coutE(InputArguments) << "RooSimPdfBuilder::buildPdf: ERROR in parsing, expected "
973  << (mode==Colon?":":"parameter list") << " in " << ruleStr->getVal() << endl ;
974  }
975 
976  //RooArgSet* paramSet = physModel->getParameters(dependents) ;
977  } else {
978  coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: no splitting rules for " << physModel->GetName() << endl ;
979  }
980  }
981 
982  coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: configured customizers for all physics models" << endl ;
983  if (oodologI((TObject*)0,ObjectHandling)) {
984  customizerList->Print() ;
985  }
986 
987  // Create fit category from physCat and splitCatList ;
988  RooArgSet fitCatList ;
989  if (physCat) fitCatList.add(*physCat) ;
990  fitCatList.add(splitCatSet) ;
991  RooSuperCategory *fitCat = new RooSuperCategory("fitCat","fitCat",fitCatList) ;
992 
993  // Create master PDF
994  RooSimultaneous* simPdf = new RooSimultaneous("simPdf","simPdf",*fitCat) ;
995 
996  // Add component PDFs to master PDF
997  TIterator* fcIter = fitCat->typeIterator() ;
998 
999  RooCatType* fcState ;
1000  while((fcState=(RooCatType*)fcIter->Next())) {
1001  // Select fitCat state
1002  fitCat->setLabel(fcState->GetName()) ;
1003 
1004  // Check if this fitCat state is selected
1005  Bool_t select(kTRUE) ;
1006  for (const auto arg : fitCatList) {
1007  auto splitCat = static_cast<const RooAbsCategory*>(arg);
1008 
1009  // Find selected state list
1010  TList* slist = (TList*) splitStateList.FindObject(splitCat->GetName()) ;
1011  if (!slist) continue ;
1012  RooCatType* type = (RooCatType*) slist->FindObject(splitCat->getCurrentLabel()) ;
1013  if (!type) {
1014  select = kFALSE ;
1015  }
1016  }
1017  if (!select) continue ;
1018 
1019 
1020  // Select appropriate PDF for this physCat state
1021  RooCustomizer* physCustomizer ;
1022  if (physCat) {
1023  RooStringVar* physNameVar = (RooStringVar*) stateMap.find(physCat->getCurrentLabel()) ;
1024  if (!physNameVar) continue ;
1025  physCustomizer = (RooCustomizer*) customizerList->FindObject(physNameVar->getVal());
1026  } else {
1027  physCustomizer = (RooCustomizer*) customizerList->First() ;
1028  }
1029 
1030  coutI(ObjectHandling) << "RooSimPdfBuilder::buildPdf: Customizing physics model " << physCustomizer->GetName()
1031  << " for mode " << fcState->GetName() << endl ;
1032 
1033  // Customizer PDF for current state and add to master simPdf
1034  RooAbsPdf* fcPdf = (RooAbsPdf*) physCustomizer->build(masterSplitCat.getCurrentLabel(),verbose) ;
1035  simPdf->addPdf(*fcPdf,fcState->GetName()) ;
1036  }
1037  delete fcIter ;
1038 
1039  // Move customizers (owning the cloned branch node components) to the attic
1040  _retiredCustomizerList.AddAll(customizerList) ;
1041  delete customizerList ;
1042 
1043  splitStateList.Delete() ;
1044 
1045  if (auxSplitCloneSet) delete auxSplitCloneSet ;
1046 
1047  _simPdfList.push_back(simPdf) ;
1048  _fitCatList.push_back(fitCat) ;
1049  return simPdf ;
1050 }
1051 
1052 
1053 
1054 
1055 
1056 ////////////////////////////////////////////////////////////////////////////////
1057 
1059 {
1061 
1062  std::list<RooSimultaneous*>::iterator iter = _simPdfList.begin() ;
1063  while(iter != _simPdfList.end()) {
1064  delete *iter ;
1065  ++iter ;
1066  }
1067 
1068  std::list<RooSuperCategory*>::iterator iter2 = _fitCatList.begin() ;
1069  while(iter2 != _fitCatList.end()) {
1070  delete *iter2 ;
1071  ++iter2 ;
1072  }
1073 
1074 }
1075 
1076 
RooFormulaVar.h
RooSuperCategory.h
RooTruthModel.h
first
Definition: first.py:1
RooHelpers.h
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:91
RooStringVar.h
RooStringVar::getVal
const char * getVal() const
Definition: RooStringVar.h:48
RooAddModel.h
TCollection::Print
virtual void Print(Option_t *option="") const
Default print for collections, calls Print(option, 1).
Definition: TCollection.cxx:476
RooMsgService.h
RooCustomizer::splitArgs
void splitArgs(const RooArgSet &argSet, const RooAbsCategory &splitCat)
Split all arguments in 'set' into individualized clones for each defined state of 'splitCat'.
Definition: RooCustomizer.cxx:293
RooFit.h
RooSimultaneous.h
RooFit::InputArguments
@ InputArguments
Definition: RooGlobalFunc.h:68
RooArgSet.h
TList::FindObject
virtual TObject * FindObject(const char *name) const
Find an object in this list using its name.
Definition: TList.cxx:577
TList::Delete
virtual void Delete(Option_t *option="")
Remove all objects from the list AND delete all heap based objects.
Definition: TList.cxx:469
RooSimPdfBuilder::RooSimPdfBuilder
RooSimPdfBuilder(const RooArgSet &pdfProtoList)
Definition: RooSimPdfBuilder.cxx:465
ClassImp
#define ClassImp(name)
Definition: Rtypes.h:364
Form
char * Form(const char *fmt,...)
RooCustomizer
Definition: RooCustomizer.h:35
RooSimPdfBuilder::_splitNodeListOwned
RooArgSet _splitNodeListOwned
Definition: RooSimPdfBuilder.h:67
coutE
#define coutE(a)
Definition: RooMsgService.h:33
coutW
#define coutW(a)
Definition: RooMsgService.h:32
RooCatType::GetName
virtual const Text_t * GetName() const
Returns name of object.
Definition: RooCatTypeLegacy.h:59
RooArgList
Definition: RooArgList.h:21
ROOT::Math::MCIntegration::Mode
Mode
Definition: MCIntegrationTypes.h:106
RooThresholdCategory.h
RooAbsCollection::remove
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE)
Remove the specified argument from our list.
Definition: RooAbsCollection.cxx:584
RooAbsCollection::find
RooAbsArg * find(const char *name) const
Find object with given name in list.
Definition: RooAbsCollection.cxx:813
RooSimPdfBuilder
Definition: RooSimPdfBuilder.h:32
coutI
#define coutI(a)
Definition: RooMsgService.h:30
RooAddPdf.h
RooArgSet::add
virtual Bool_t add(const RooAbsCollection &col, Bool_t silent=kFALSE)
Add a collection of arguments to this collection by calling add() for each element in the source coll...
Definition: RooArgSet.h:88
RooAbsReal
Definition: RooAbsReal.h:61
RooAbsArg::recursiveRedirectServers
Bool_t recursiveRedirectServers(const RooAbsCollection &newServerList, Bool_t mustReplaceAll=kFALSE, Bool_t nameChange=kFALSE, Bool_t recurseInNewSet=kTRUE)
Definition: RooAbsArg.cxx:1066
TString
Definition: TString.h:136
TCollection::AddAll
virtual void AddAll(const TCollection *col)
Add all objects from collection col to this collection.
Definition: TCollection.cxx:195
RooGenericPdf.h
RooSimPdfBuilder::createProtoBuildConfig
RooArgSet * createProtoBuildConfig()
Make RooArgSet of configuration objects.
Definition: RooSimPdfBuilder.cxx:479
RooCustomizer::build
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...
Definition: RooCustomizer.cxx:395
RooAbsCategory::getCurrentLabel
virtual const char * getCurrentLabel() const
Return label string of current state.
Definition: RooAbsCategory.cxx:130
RooSimultaneous::addPdf
Bool_t addPdf(const RooAbsPdf &pdf, const char *catLabel)
Associate given PDF with index category state label 'catLabel'.
Definition: RooSimultaneous.cxx:373
RooDataSet.h
RooSimPdfBuilder::_retiredCustomizerList
TList _retiredCustomizerList
Definition: RooSimPdfBuilder.h:69
RooAbsCollection::useHashMapForFind
void useHashMapForFind(bool flag) const
Install an internal hash map for fast finding of elements by name.
Definition: RooAbsCollection.cxx:1314
bool
TIterator
Definition: TIterator.h:30
RooArgSet::addOwned
virtual Bool_t addOwned(const RooAbsCollection &col, Bool_t silent=kFALSE)
Add a collection of arguments to this collection by calling addOwned() for each element in the source...
Definition: RooArgSet.h:92
RooAbsCategory
Definition: RooAbsCategory.h:38
RooAbsArg::getParameters
RooArgSet * getParameters(const RooAbsData *data, Bool_t stripDisconnected=kTRUE) 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:546
RooArgSet::snapshot
RooArgSet * snapshot(bool deepCopy=true) const
Use RooAbsCollection::snapshot(), but return as RooArgSet.
Definition: RooArgSet.h:134
TList::First
virtual TObject * First() const
Return the first object in the list. Returns 0 when list is empty.
Definition: TList.cxx:658
TObject::GetName
virtual const char * GetName() const
Returns name of object.
Definition: TObject.cxx:359
RooAbsArg::clone
virtual TObject * clone(const char *newname=0) const =0
RooTrace.h
RooLinearVar.h
RooFormulaVar
Definition: RooFormulaVar.h:30
RooCustomizer.h
RooSimPdfBuilder::_fitCatList
std::list< RooSuperCategory * > _fitCatList
Definition: RooSimPdfBuilder.h:72
RooProdPdf.h
RooSimPdfBuilder::addSpecializations
void addSpecializations(const RooArgSet &specSet)
Definition: RooSimPdfBuilder.cxx:499
RooMultiCategory
Definition: RooMultiCategory.h:28
RooMultiCategory.h
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:92
TString::Append
TString & Append(const char *cs)
Definition: TString.h:564
RooSimPdfBuilder::buildPdf
RooSimultaneous * buildPdf(const RooArgSet &buildConfig, const RooArgSet &dependents, const RooArgSet *auxSplitCats=0, Bool_t verbose=kFALSE)
Initialize needed components.
Definition: RooSimPdfBuilder.cxx:509
RooDataHist.h
RooPlot.h
RooAbsCollection::createIterator
TIterator * createIterator(Bool_t dir=kIterForward) const
TIterator-style iteration over contained elements.
Definition: RooAbsCollection.h:118
oodologI
#define oodologI(o, a)
Definition: RooMsgService.h:73
RooAbsCollection::add
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
Definition: RooAbsCollection.cxx:437
RooSimPdfBuilder::_protoPdfSet
RooArgSet _protoPdfSet
Definition: RooSimPdfBuilder.h:64
RooAbsCategory::lookupType
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.
Definition: RooAbsCategory.cxx:307
RooCategory.h
RooAbsCategory::typeIterator
TIterator * typeIterator() const
Definition: RooAbsCategory.cxx:653
RooFit::ObjectHandling
@ ObjectHandling
Definition: RooGlobalFunc.h:68
RooSuperCategory::setLabel
virtual Bool_t setLabel(const char *label, Bool_t printError=kTRUE) override
Set the value of the super category by specifying the state name.
Definition: RooSuperCategory.cxx:132
RooRealVar.h
RooSimPdfBuilder::_simPdfList
std::list< RooSimultaneous * > _simPdfList
Definition: RooSimPdfBuilder.h:71
RooFitResult.h
TIterator::Next
virtual TObject * Next()=0
RooSimPdfBuilder::_compSplitCatSet
RooArgSet _compSplitCatSet
Definition: RooSimPdfBuilder.h:66
RooHelpers::tokenise
std::vector< std::string > tokenise(const std::string &str, const std::string &delims, bool returnEmptyToken=true)
Tokenise the string by splitting at the characters in delims.
Definition: RooHelpers.cxx:74
proto
const char * proto
Definition: civetweb.c:16604
RooSimPdfBuilder::~RooSimPdfBuilder
~RooSimPdfBuilder()
Definition: RooSimPdfBuilder.cxx:1058
TString::CompareTo
int CompareTo(const char *cs, ECaseCompare cmp=kExact) const
Compare a string to char *cs2.
Definition: TString.cxx:418
RooAbsCategory.h
RooSimPdfBuilder::_splitNodeList
RooArgSet _splitNodeList
Definition: RooSimPdfBuilder.h:68
RooStringVar
Definition: RooStringVar.h:23
RooCategory
Definition: RooCategory.h:27
TCollection::SetName
void SetName(const char *name)
Definition: TCollection.h:204
RooSuperCategory
Definition: RooSuperCategory.h:27
TList::Add
virtual void Add(TObject *obj)
Definition: TList.h:87
TObject
Definition: TObject.h:37
RooAbsCategory::hasLabel
bool hasLabel(const std::string &label) const
Check if a state with name label exists.
Definition: RooAbsCategory.h:69
RooCatType
Definition: RooCatTypeLegacy.h:23
genreflex::verbose
bool verbose
Definition: rootcling_impl.cxx:133
RooAbsArg
Definition: RooAbsArg.h:73
RooRealIntegral.h
RooAbsPdf
Definition: RooAbsPdf.h:40
TNamed::GetName
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:53
type
int type
Definition: TGX11.cxx:121
RooSimultaneous
Definition: RooSimultaneous.h:37
RooAbsCategoryLValue
Definition: RooAbsCategoryLValue.h:25
Riostream.h
RooSimPdfBuilder.h
RooMappedCategory.h
RooAbsCollection::getSize
Int_t getSize() const
Definition: RooAbsCollection.h:171
TList
Definition: TList.h:44
RooArgSet
Definition: RooArgSet.h:28
int