Logo ROOT  
Reference Guide
 
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
Loading...
Searching...
No Matches
SPlot.cxx
Go to the documentation of this file.
1// @(#)root/roostats:$Id$
2// Author: Kyle Cranmer 28/07/2008
3
4/*************************************************************************
5 * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12/*****************************************************************************
13 * Project: RooStats
14 * Package: RooFit/RooStats
15 *
16 * Authors:
17 * Original code from M. Pivk as part of MLFit package from BaBar.
18 * Modifications:
19 * Giacinto Piacquadio, Maurizio Pierini: modifications for new RooFit version
20 * George H. Lewis, Kyle Cranmer: generalized for weighted events
21 *
22 * Porting to RooStats (with permission) by Kyle Cranmer, July 2008
23 *
24 *****************************************************************************/
25
26
27/** \class RooStats::SPlot
28 \ingroup Roostats
29
30 A class to calculate "sWeights" used to create an "sPlot".
31 An sPlot can reweight a dataset to show different components (e.g. signal / background),
32 but it doesn't use cuts, and therefore doesn't have to sort events into signal/background (or other) categories.
33 Instead of *assigning* a category to each event in the dataset, all events are *weighted*.
34 To compute the weights, a PDF with different components is analysed, and the weights are added
35 to the dataset. When plotting the dataset with the weights of the signal or
36 background components, the data looks like "signal", but all events in the dataset are used.
37
38 The result is similar to a likelihood projection plot, but without cuts.
39
40 \note SPlot needs to fit the pdf to the data once, so make sure that all relevant fit arguments such as
41 the fit range are passed in the constructor.
42
43 The code is based on
44 ``SPlot: A statistical tool to unfold data distributions,''
45 Nucl. Instrum. Meth. A 555, 356 (2005) [arXiv:physics/0402083].
46
47 ### Creating an SPlot
48 To use this class, you first must have a pdf that includes
49 yield parameters for (possibly several) different species, for example a signal and background
50 yield. Those yields must be of type RooRealVar / RooLinearVar (or anything that derives from
51 RooAbsRealLValue). This is necessary because
52 RooStats needs to be able to set the yields to 0 and 1 to probe the PDF. After
53 constructing the s weights, the yields will be restored to their original values.
54
55 To create an instance of the SPlot, supply a data set, the pdf to analyse,
56 and a list which parameters of the pdf are yields. The SPlot will calculate SWeights, and
57 include these as columns in the RooDataSet. The dataset will have two additional columns
58 for every yield with name "`<varname>`":
59 - `L_<varname>` is the likelihood for each event, *i.e.*, the pdf evaluated for the given value of the variable "varname".
60 - `<varname>_sw` is the value of the sWeight for the variable "varname" for each event.
61
62 In SPlot::SPlot(), one can choose whether columns should be added to an existing dataset or whether a copy of the dataset
63 should be created.
64
65 ### Plotting s-weighted data
66 After computing the s weights, create a new dataset that uses the s weights of the variable of interest for weighting.
67 If the yield parameter for signal was e.g. "signalYield", the dataset can be constructed as follows:
68 ~~~{.cpp}
69 RooDataSet data_signal("<name>", "<title>", <dataWithSWeights>, <variables>, 0, "signalYield_sw");
70 ~~~
71
72 A complete tutorial with an extensive model is rs301_splot.C
73
74 #### Using ratios as yield parameters
75 As mentioned, RooStats needs to be able to modify the yield parameters. That means that they have to be a RooRealVar
76 of a RooLinearVar. This allows using ratio parameters as in the following example:
77 ~~~{.cpp}
78 RooRealVar x("x", "observable", 0, 0, 20);
79 RooRealVar m("m", "mean", 5., -10, 10);
80 RooRealVar s("s", "sigma", 2., 0, 10);
81 RooGaussian gauss("gauss", "gauss", x, m, s);
82
83 RooRealVar a("a", "exp", -0.2, -10., 0.);
84 RooExponential ex("ex", "ex", x, a);
85
86 RooRealVar common("common", "common scale", 3., 0, 10);
87 RooRealVar r1("r1", "ratio of signal events", 0.3, 0, 10);
88 RooRealVar r2("r2", "ratio of background events", 0.5, 0, 10);
89 RooLinearVar c1("c1", "c1", r1, common, RooFit::RooConst(0.));
90 RooLinearVar c2("c2", "c2", r2, common, RooFit::RooConst(0.));
91
92 RooAddPdf sum("sum", "sum", RooArgSet(gauss, ex), RooArgSet(c1, c2));
93 auto data = sum.generate(x, 1000);
94
95 RooStats::SPlot splot("splot", "splot", *data, &sum, RooArgSet(c1, c2));
96 ~~~
97*/
98
99#include <vector>
100#include <map>
101
102#include "RooStats/SPlot.h"
103#include "RooAbsPdf.h"
104#include "RooDataSet.h"
105#include "RooRealVar.h"
106#include "RooGlobalFunc.h"
108
109
110#include "TMatrixD.h"
111
112
113
114using namespace RooStats;
115using std::endl;
116
117////////////////////////////////////////////////////////////////////////////////
118
120{
121 if(TestBit(kOwnData) && fSData)
122 delete fSData;
123
124}
125
126////////////////////////////////////////////////////////////////////////////////
127/// Default constructor
128
130{
131 RooArgList Args;
132
133 fSWeightVars.assign(Args);
134}
135
136////////////////////////////////////////////////////////////////////////////////
137
138SPlot::SPlot(const char *name, const char *title) : TNamed(name, title)
139{
140 RooArgList Args;
141 fSWeightVars.assign(Args);
142}
143
144////////////////////////////////////////////////////////////////////////////////
145///Constructor from a RooDataSet
146///No sWeighted variables are present
147
148SPlot::SPlot(const char *name, const char *title, const RooDataSet &data)
149 : TNamed(name, title), fSData(const_cast<RooDataSet *>(&data))
150{
151 RooArgList Args;
152
153 fSWeightVars.assign(Args);
154}
155
156////////////////////////////////////////////////////////////////////////////////
157/// Copy Constructor from another SPlot
158
159SPlot::SPlot(const SPlot &other) : TNamed(other), fSData((RooDataSet *)other.GetSDataSet())
160{
161 RooArgList Args = (RooArgList) other.GetSWeightVars();
162
164}
165
166////////////////////////////////////////////////////////////////////////////////
167///Construct a new SPlot instance, calculate sWeights, and include them
168///in the RooDataSet held by this instance.
169///
170/// The constructor automatically calls AddSWeight() to add s weights to the dataset.
171/// These can be retrieved later using GetSWeight() or GetSDataSet().
172///\param[in] name Name of the instance.
173///\param[in] title Title of the instance.
174///\param[in] data Dataset to fit to.
175///\param[in] pdf PDF to compute s weights for.
176///\param[in] yieldsList List of parameters in `pdf` that are yields. These must be RooRealVar or RooLinearVar, since RooStats will need to modify their values.
177///\param[in] projDeps Don't normalise over these parameters when calculating the sWeights. Will be passed on to AddSWeight().
178///\param[in] useWeights Include weights of the input data in calculation of s weights.
179///\param[in] cloneData Make a clone of the incoming data before adding weights.
180///\param[in] newName New name for the data.
181///\param[in] arg5,arg6,arg7,arg8 Additional arguments for the fitting step in AddSWeight().
182SPlot::SPlot(const char* name, const char* title, RooDataSet& data, RooAbsPdf* pdf,
183 const RooArgList &yieldsList, const RooArgSet &projDeps,
184 bool useWeights, bool cloneData, const char* newName,
185 const RooCmdArg& arg5, const RooCmdArg& arg6, const RooCmdArg& arg7, const RooCmdArg& arg8):
186 TNamed(name, title)
187{
188 if(cloneData) {
189 fSData = static_cast<RooDataSet*>(data.Clone(newName));
191 }
192 else
193 fSData = (RooDataSet*) &data;
194
195 // Add check that yieldsList contains all RooRealVar / RooAbsRealLValue
196 for (const auto arg : yieldsList) {
197 if (!dynamic_cast<const RooAbsRealLValue*>(arg)) {
198 coutE(InputArguments) << "SPlot::SPlot(" << GetName() << ") input argument "
199 << arg->GetName() << " is not of type RooRealVar (or RooLinearVar)."
200 << "\nRooStats must be able to set it to 0 and to 1 to probe the PDF." << std::endl ;
201 throw std::invalid_argument(Form("SPlot::SPlot(%s) input argument %s is not of type RooRealVar/RooLinearVar",GetName(),arg->GetName())) ;
202 }
203 }
204
205 //Construct a new SPlot class,
206 //calculate sWeights, and include them
207 //in the RooDataSet of this class.
208
209 this->AddSWeight(pdf, yieldsList, projDeps, useWeights, arg5, arg6, arg7, arg8);
210}
211
212////////////////////////////////////////////////////////////////////////////////
213/// Set dataset (if not passed in constructor).
215{
216 if(data) {
218 return fSData;
219 } else
220 return nullptr;
221}
222
223////////////////////////////////////////////////////////////////////////////////
224/// Retrieve s-weighted data.
225/// It does **not** automatically call AddSWeight(). This needs to be done manually.
227{
228 return fSData;
229}
230
231////////////////////////////////////////////////////////////////////////////////
232/// Retrieve an s weight.
233/// \param[in] numEvent Event number to retrieve s weight for.
234/// \param[in] sVariable The yield parameter to retrieve the s weight for.
235double SPlot::GetSWeight(Int_t numEvent, const char* sVariable) const
236{
237 if(numEvent > fSData->numEntries() )
238 {
239 coutE(InputArguments) << "Invalid Entry Number" << std::endl;
240 return -1;
241 }
242
243 if(numEvent < 0)
244 {
245 coutE(InputArguments) << "Invalid Entry Number" << std::endl;
246 return -1;
247 }
248
249 double totalYield = 0;
250
251 std::string varname(sVariable);
252 varname += "_sw";
253
254
256 {
259
260 return totalYield;
261 }
262
263 if( fSWeightVars.find(varname.c_str()) )
264 {
265
267 totalYield += Row.getRealValue(varname.c_str() );
268
269 return totalYield;
270 }
271
272 else
273 coutE(InputArguments) << "InputVariable not in list of sWeighted variables" << std::endl;
274
275 return -1;
276}
277
278
279////////////////////////////////////////////////////////////////////////////////
280/// Sum the SWeights for a particular event.
281/// This sum should equal the total weight of that event.
282/// This method is intended to be used as a check.
283
285{
286 if(numEvent > fSData->numEntries() )
287 {
288 coutE(InputArguments) << "Invalid Entry Number" << std::endl;
289 return -1;
290 }
291
292 if(numEvent < 0)
293 {
294 coutE(InputArguments) << "Invalid Entry Number" << std::endl;
295 return -1;
296 }
297
299
300 double eventSWeight = 0;
301
303
304 for (Int_t i = 0; i < numSWeightVars; i++)
306
307 return eventSWeight;
308}
309
310////////////////////////////////////////////////////////////////////////////////
311/// Sum the SWeights for a particular species over all events.
312/// This should equal the total (weighted) yield of that species.
313/// This method is intended as a check.
314
315double SPlot::GetYieldFromSWeight(const char* sVariable) const
316{
317
318 double totalYield = 0;
319
320 std::string varname(sVariable);
321 varname += "_sw";
322
323
325 {
326 for(Int_t i=0; i < fSData->numEntries(); i++)
327 {
328 RooArgSet Row(*fSData->get(i));
330 }
331
332 return totalYield;
333 }
334
335 if( fSWeightVars.find(varname.c_str()) )
336 {
337 for(Int_t i=0; i < fSData->numEntries(); i++)
338 {
339 RooArgSet Row(*fSData->get(i));
340 totalYield += Row.getRealValue(varname.c_str() );
341 }
342
343 return totalYield;
344 }
345
346 else
347 coutE(InputArguments) << "InputVariable not in list of sWeighted variables" << std::endl;
348
349 return -1;
350}
351
352
353////////////////////////////////////////////////////////////////////////////////
354/// Return a RooArgList containing all parameters that have s weights.
355
357{
358
360
361 return Args;
362
363}
364
365////////////////////////////////////////////////////////////////////////////////
366/// Return the number of SWeights
367/// In other words, return the number of
368/// species that we are trying to extract.
369
371{
373
374 return Args.size();
375}
376
377////////////////////////////////////////////////////////////////////////////////
378/// Method which adds the sWeights to the dataset.
379///
380/// The SPlot will contain two new variables for each yield parameter:
381/// - `L_<varname>` is the likelihood for each event, i.e., the pdf evaluated for the a given value of the variable "varname".
382/// - `<varname>_sw` is the value of the sWeight for the variable "varname" for each event.
383///
384/// Find Parameters in the PDF to be considered fixed when calculating the SWeights
385/// and be sure to NOT include the yields in that list.
386///
387/// After fixing non-yield parameters, this function will start a fit by calling
388/// ```
389/// pdf->fitTo(*fSData, RooFit::Extended(true), RooFit::SumW2Error(true), RooFit::PrintLevel(-1), RooFit::PrintEvalErrors(-1)).
390/// ```
391/// One can pass additional arguments to `fitTo`, such as `RooFit::Range("fitrange")`, as `arg5`, `arg6`, `arg7`, `arg8`.
392///
393/// \note A `RooFit::Range` may be necessary to get expected results if you initially fit in a range
394/// and/or called `pdf->fixCoefRange("fitrange")` on `pdf`.
395/// Pass `arg5`, `arg6`, `arg7`, `arg8` AT YOUR OWN RISK.
396///
397/// \param[in] pdf PDF to fit to data to compute s weights.
398/// \param[in] yieldsTmp Yields to use to compute s weights.
399/// \param[in] projDeps These will not be normalized over when calculating the sWeights,
400/// and will be considered parameters, not observables.
401/// \param[in] includeWeights Include weights of the input data in calculation of s weights.
402/// \param[in] arg5,arg6,arg7,arg8 Optional additional arguments for the fitting step.
404 const RooArgSet &projDeps, bool includeWeights,
405 const RooCmdArg& arg5, const RooCmdArg& arg6, const RooCmdArg& arg7, const RooCmdArg& arg8)
406{
407
408 // Find Parameters in the PDF to be considered fixed when calculating the SWeights
409 // and be sure to NOT include the yields in that list
410 std::unique_ptr<RooArgSet> constParameters{pdf->getParameters(fSData)};
411 for (unsigned int i=0; i < constParameters->size(); ++i) {
412 // Need a counting loop since collection is being modified
413 auto& par = *(*constParameters)[i];
414 if (std::any_of(yieldsTmp.begin(), yieldsTmp.end(), [&](const RooAbsArg* yield){ return yield->dependsOn(par); })) {
415 constParameters->remove(par, true, true);
416 --i;
417 }
418 }
419
420
421 // Set these parameters constant and store them so they can later
422 // be set to not constant
423 std::vector<RooAbsRealLValue*> constVarHolder;
424
425 for(std::size_t i = 0; i < constParameters->size(); i++)
426 {
427 RooAbsRealLValue* varTemp = static_cast<RooAbsRealLValue*>((*constParameters)[i] );
428 if(varTemp && varTemp->isConstant() == 0 )
429 {
430 varTemp->setConstant();
431 constVarHolder.push_back(varTemp);
432 }
433 }
434
435 // Fit yields to the data with all other variables held constant
436 // This is necessary because SPlot assumes the yields minimise -Log(likelihood)
438
439 // The list of variables to normalize over when calculating PDF values.
440 RooArgSet vars(*fSData->get() );
441 vars.remove(projDeps, true, true);
442
443 // Hold the value of the fitted yields
444 std::vector<double> yieldsHolder;
445
446 yieldsHolder.reserve(yieldsTmp.size());
447 for(std::size_t i = 0; i < yieldsTmp.size(); i++) {
448 yieldsHolder.push_back(static_cast<RooAbsReal*>(yieldsTmp.at(i))->getVal(&vars));
449 }
450
451 const Int_t nspec = yieldsTmp.size();
452 RooArgList yields = *static_cast<RooArgList*>(yieldsTmp.snapshot(false));
453
455 coutI(InputArguments) << "Printing Yields" << std::endl;
456 yields.Print();
457 }
458
459
460 // first calculate the pdf values for all species and all events
461 std::vector<RooAbsRealLValue*> yieldvars ;
464
465 std::vector<double> yieldvalues ;
466 for (Int_t k = 0; k < nspec; ++k) {
467 auto thisyield = static_cast<const RooAbsReal*>(yields.at(k)) ;
468 auto yieldinpdf = static_cast<RooAbsRealLValue*>(pdfServers.find(thisyield->GetName()));
470
471 if (yieldinpdf) {
472 coutI(InputArguments)<< "yield in pdf: " << yieldinpdf->GetName() << " " << thisyield->getVal(&vars) << std::endl;
473
474 yieldvars.push_back(yieldinpdf) ;
475 yieldvalues.push_back(thisyield->getVal(&vars)) ;
476 }
477 }
478
480
481
482
483
484 // set all yield to zero
485 for(Int_t m=0; m<nspec; ++m) {
486 auto theVar = static_cast<RooAbsRealLValue*>(yieldvars[m]);
487 theVar->setVal(0) ;
488
489 //Check that range of yields is at least (0,1), and fix otherwise
490 if (theVar->getMin() > 0) {
491 coutE(InputArguments) << "Yield variables need to have a range that includes at least [0, 1]. Minimum for "
492 << theVar->GetName() << " is " << theVar->getMin() << std::endl;
493 if (RooRealVar* realVar = dynamic_cast<RooRealVar*>(theVar)) {
494 coutE(InputArguments) << "Setting min range to 0" << std::endl;
495 realVar->setMin(0);
496 } else {
497 throw std::invalid_argument(std::string("Yield variable ") + theVar->GetName() + " must have a range that includes 0.");
498 }
499 }
500
501 if (theVar->getMax() < 1) {
502 coutW(InputArguments) << "Yield variables need to have a range that includes at least [0, 1]. Maximum for "
503 << theVar->GetName() << " is " << theVar->getMax() << std::endl;
504 if (RooRealVar* realVar = dynamic_cast<RooRealVar*>(theVar)) {
505 coutE(InputArguments) << "Setting max range to 1" << std::endl;
506 realVar->setMax(1);
507 } else {
508 throw std::invalid_argument(std::string("Yield variable ") + theVar->GetName() + " must have a range that includes 1.");
509 }
510 }
511 }
512
513
514 // For every event and for every species,
515 // calculate the value of the component pdf for that specie
516 // by setting the yield of that specie to 1
517 // and all others to 0. Evaluate the pdf for each event
518 // and store the values.
519
520 std::unique_ptr<RooArgSet> pdfvars{pdf->getVariables()};
521 std::vector<std::vector<double> > pdfvalues(numevents,std::vector<double>(nspec,0)) ;
522
523 for (Int_t ievt = 0; ievt <numevents; ievt++)
524 {
525 //WVE: FIX THIS PART, EVALUATION PROGRESS!!
526
527 pdfvars->assign(*fSData->get(ievt));
528
529 for(Int_t k = 0; k < nspec; ++k) {
530 auto theVar = static_cast<RooAbsRealLValue*>(yieldvars[k]);
531
532 // set this yield to 1
533 theVar->setVal( 1 ) ;
534 // evaluate the pdf
535 double f_k = pdf->getVal(&vars) ;
536 pdfvalues[ievt][k] = f_k ;
537 if( !(f_k>1 || f_k<1) )
538 coutW(InputArguments) << "Strange pdf value: " << ievt << " " << k << " " << f_k << std::endl ;
539 theVar->setVal( 0 ) ;
540 }
541 }
542
543 // check that the likelihood normalization is fine
544 std::vector<double> norm(nspec,0) ;
545 for (Int_t ievt = 0; ievt <numevents ; ievt++)
546 {
547 double dnorm(0) ;
548 for(Int_t k=0; k<nspec; ++k) dnorm += yieldvalues[k] * pdfvalues[ievt][k] ;
549 for(Int_t j=0; j<nspec; ++j) norm[j] += pdfvalues[ievt][j]/dnorm ;
550 }
551
552 coutI(Contents) << "likelihood norms: " ;
553
554 for(Int_t k=0; k<nspec; ++k) coutI(Contents) << norm[k] << " " ;
555 coutI(Contents) << std::endl ;
556
557 // Make a TMatrixD to hold the covariance matrix.
559 for (Int_t i = 0; i < nspec; i++) for (Int_t j = 0; j < nspec; j++) covInv(i,j) = 0;
560
561 coutI(Contents) << "Calculating covariance matrix";
562
563
564 // Calculate the inverse covariance matrix, using weights
565 for (Int_t ievt = 0; ievt < numevents; ++ievt)
566 {
567
568 fSData->get(ievt) ;
569
570 // Calculate contribution to the inverse of the covariance
571 // matrix. See BAD 509 V2 eqn. 15
572
573 // Sum for the denominator
574 double dsum(0);
575 for(Int_t k = 0; k < nspec; ++k)
576 dsum += pdfvalues[ievt][k] * yieldvalues[k] ;
577
578 for (Int_t n = 0; n < nspec; ++n) {
579 for(Int_t j=0; j<nspec; ++j)
580 {
581 if (includeWeights) {
582 covInv(n, j) += fSData->weight() * pdfvalues[ievt][n] * pdfvalues[ievt][j] / (dsum * dsum);
583 } else {
584 covInv(n, j) += pdfvalues[ievt][n] * pdfvalues[ievt][j] / (dsum * dsum);
585 }
586 }
587 }
588
589 //ADDED WEIGHT ABOVE
590
591 }
592
593 // Covariance inverse should now be computed!
594
595 // Invert to get the covariance matrix
596 if (covInv.Determinant() <=0)
597 {
598 coutE(Eval) << "SPlot Error: covariance matrix is singular; I can't invert it!" << std::endl;
599 covInv.Print();
600 return;
601 }
602
604
605 //check cov normalization
606 if (RooMsgService::instance().isActive(this, RooFit::Eval, RooFit::DEBUG)) {
607 coutI(Eval) << "Checking Likelihood normalization: " << std::endl;
608 coutI(Eval) << "Yield of specie Sum of Row in Matrix Norm" << std::endl;
609 for(Int_t k=0; k<nspec; ++k)
610 {
611 double covnorm(0) ;
612 for(Int_t m=0; m<nspec; ++m) covnorm += covInv[k][m]*yieldvalues[m] ;
613 double sumrow(0) ;
614 for(Int_t m = 0; m < nspec; ++m) sumrow += covMatrix[k][m] ;
615 coutI(Eval) << yieldvalues[k] << " " << sumrow << " " << covnorm << std::endl ;
616 }
617 }
618
619 // calculate for each event the sWeight (BAD 509 V2 eq. 21)
620 coutI(Eval) << "Calculating sWeight" << std::endl;
621 std::vector<RooRealVar*> sweightvec ;
622 std::vector<RooRealVar*> pdfvec ;
624
625 // Create and label the variables
626 // used to store the SWeights
627
629
630 for(Int_t k=0; k<nspec; ++k)
631 {
632 std::string wname = std::string(yieldvars[k]->GetName()) + "_sw";
633 RooRealVar* var = new RooRealVar(wname.c_str(),wname.c_str(),0) ;
634 sweightvec.push_back( var) ;
635 sweightset.add(*var) ;
636 fSWeightVars.add(*var);
637
638 wname = "L_" + std::string(yieldvars[k]->GetName());
639 var = new RooRealVar(wname.c_str(),wname.c_str(),0) ;
640 pdfvec.push_back( var) ;
641 sweightset.add(*var) ;
642 }
643
644 // Create and fill a RooDataSet
645 // with the SWeights
646
647 RooDataSet* sWeightData = new RooDataSet("dataset", "dataset with sWeights", sweightset);
648
649 for(Int_t ievt = 0; ievt < numevents; ++ievt)
650 {
651
652 fSData->get(ievt) ;
653
654 // sum for denominator
655 double dsum(0);
656 for(Int_t k = 0; k < nspec; ++k) dsum += pdfvalues[ievt][k] * yieldvalues[k] ;
657 // covariance weighted pdf for each specief
658 for(Int_t n=0; n<nspec; ++n)
659 {
660 double nsum(0) ;
661 for(Int_t j=0; j<nspec; ++j) nsum += covMatrix(n,j) * pdfvalues[ievt][j] ;
662
663
664 //Add the sWeights here!!
665 //Include weights,
666 //ie events weights are absorbed into sWeight
667
668
669 if(includeWeights) sweightvec[n]->setVal(fSData->weight() * nsum/dsum) ;
670 else sweightvec[n]->setVal( nsum/dsum) ;
671
672 pdfvec[n]->setVal( pdfvalues[ievt][n] ) ;
673
674 if( !(std::abs(nsum/dsum)>=0 ) )
675 {
676 coutE(Contents) << "error: " << nsum/dsum << std::endl ;
677 return;
678 }
679 }
680
681 sWeightData->add(sweightset) ;
682 }
683
684
685 // Add the SWeights to the original data set
686
688
689 delete sWeightData;
690
691 //Restore yield values
692
693 for(std::size_t i = 0; i < yieldsTmp.size(); i++)
694 static_cast<RooAbsRealLValue*>(yieldsTmp.at(i))->setVal(yieldsHolder.at(i));
695
696 //Make any variables that were forced to constant no longer constant
697
698 for(Int_t i=0; i < (Int_t) constVarHolder.size(); i++)
699 constVarHolder.at(i)->setConstant(false);
700
701 return;
702
703}
#define coutI(a)
#define coutW(a)
#define coutE(a)
int Int_t
Definition RtypesCore.h:45
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
char name[80]
Definition TGX11.cxx:110
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition TString.cxx:2489
const_iterator begin() const
const_iterator end() const
Common abstract base class for objects that represent a value and a "shape" in RooFit.
Definition RooAbsArg.h:77
bool dependsOn(const RooAbsCollection &serverList, const RooAbsArg *ignoreArg=nullptr, bool valueOnly=false) const
Test whether we depend on (ie, are served by) any object in the specified collection.
RooFit::OwningPtr< 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...
RooFit::OwningPtr< RooArgSet > getVariables(bool stripDisconnected=true) const
Return RooArgSet with all variables (tree leaf nodes of expression tree)
void treeNodeServerList(RooAbsCollection *list, const RooAbsArg *arg=nullptr, bool doBranch=true, bool doLeaf=true, bool valueOnly=false, bool recurseNonDerived=false) const
Fill supplied list with nodes of the arg tree, following all server links, starting with ourself as t...
double getRealValue(const char *name, double defVal=0.0, bool verbose=false) const
Get value of a RooAbsReal stored in set with given name.
virtual bool remove(const RooAbsArg &var, bool silent=false, bool matchByNameOnly=false)
Remove the specified argument from our list.
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
void assign(const RooAbsCollection &other) const
Sets the value, cache and constant attribute of any argument in our set that also appears in the othe...
Storage_t::size_type size() const
virtual RooAbsArg * addClone(const RooAbsArg &var, bool silent=false)
Add a clone of the specified argument to list.
RooAbsArg * find(const char *name) const
Find object with given name in list.
virtual Int_t numEntries() const
Return number of entries in dataset, i.e., count unweighted entries.
Abstract interface for all probability density functions.
Definition RooAbsPdf.h:40
RooFit::OwningPtr< RooFitResult > fitTo(RooAbsData &data, CmdArgs_t const &... cmdArgs)
Fit PDF to given dataset.
Definition RooAbsPdf.h:157
Abstract base class for objects that represent a real value that may appear on the left hand side of ...
virtual void setVal(double value)=0
Set the current value of the object. Needs to be overridden by implementations.
Abstract base class for objects that represent a real value and implements functionality common to al...
Definition RooAbsReal.h:59
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition RooAbsReal.h:103
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition RooArgList.h:22
RooAbsArg * at(Int_t idx) const
Return object at given index, or nullptr if index is out of range.
Definition RooArgList.h:110
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:24
Named container for two doubles, two integers two object points and three string pointers that can be...
Definition RooCmdArg.h:26
Container class to hold unbinned data.
Definition RooDataSet.h:34
const RooArgSet * get(Int_t index) const override
Return RooArgSet with coordinates of event 'index'.
bool merge(RooDataSet *data1, RooDataSet *data2=nullptr, RooDataSet *data3=nullptr, RooDataSet *data4=nullptr, RooDataSet *data5=nullptr, RooDataSet *data6=nullptr)
double weight() const override
Return event weight of current event.
static RooMsgService & instance()
Return reference to singleton instance.
Variable that can be changed from the outside.
Definition RooRealVar.h:37
A class to calculate "sWeights" used to create an "sPlot".
Definition SPlot.h:32
double GetSWeight(Int_t numEvent, const char *sVariable) const
Retrieve an s weight.
Definition SPlot.cxx:235
void AddSWeight(RooAbsPdf *pdf, const RooArgList &yieldsTmp, const RooArgSet &projDeps=RooArgSet(), bool includeWeights=true, const RooCmdArg &fitToarg5={}, const RooCmdArg &fitToarg6={}, const RooCmdArg &fitToarg7={}, const RooCmdArg &fitToarg8={})
Method which adds the sWeights to the dataset.
Definition SPlot.cxx:403
RooDataSet * fSData
Definition SPlot.h:82
RooArgList GetSWeightVars() const
Return a RooArgList containing all parameters that have s weights.
Definition SPlot.cxx:356
double GetSumOfEventSWeight(Int_t numEvent) const
Sum the SWeights for a particular event.
Definition SPlot.cxx:284
SPlot()
Default constructor.
Definition SPlot.cxx:129
Int_t GetNumSWeightVars() const
Return the number of SWeights In other words, return the number of species that we are trying to extr...
Definition SPlot.cxx:370
~SPlot() override
Definition SPlot.cxx:119
RooDataSet * SetSData(RooDataSet *data)
Set dataset (if not passed in constructor).
Definition SPlot.cxx:214
RooDataSet * GetSDataSet() const
Retrieve s-weighted data.
Definition SPlot.cxx:226
double GetYieldFromSWeight(const char *sVariable) const
Sum the SWeights for a particular species over all events.
Definition SPlot.cxx:315
RooArgList fSWeightVars
Definition SPlot.h:78
The TNamed class is the base class for all named ROOT classes.
Definition TNamed.h:29
const char * GetName() const override
Returns name of object.
Definition TNamed.h:47
virtual void Clear(Option_t *="")
Definition TObject.h:119
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition TObject.h:199
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition TObject.cxx:798
RooCmdArg SumW2Error(bool flag)
RooCmdArg PrintEvalErrors(Int_t numErrors)
RooCmdArg PrintLevel(Int_t code)
RooCmdArg Extended(bool flag=true)
const Int_t n
Definition legend1.C:16
@ InputArguments
Namespace for the RooStats classes.
Definition CodegenImpl.h:58
TMarker m
Definition textangle.C:8