Logo ROOT   6.16/01
Reference Guide
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 This class calculates sWeights used to create an sPlot.
31 The code is based on
32 ``SPlot: A statistical tool to unfold data distributions,''
33 Nucl. Instrum. Meth. A 555, 356 (2005) [arXiv:physics/0402083].
34
35 An SPlot gives us the distribution of some variable, x in our
36 data sample for a given species (eg. signal or background).
37 The result is similar to a likelihood projection plot, but no cuts are made,
38 so every event contributes to the distribution.
39
40 To use this class, you first must have a pdf that includes
41 yields for (possibly several) different species.
42 Create an instance of the class by supplying a data set,
43 the pdf, and a list of the yield variables. The SPlot Class
44 will calculate SWeights and include these as columns in the RooDataSet.
45
46*/
47
48#include <vector>
49#include <map>
50
51#include "RooStats/SPlot.h"
52#include "RooAbsPdf.h"
53#include "RooDataSet.h"
54#include "RooRealVar.h"
55#include "RooGlobalFunc.h"
56#include "TTree.h"
58
59
60#include "TMatrixD.h"
61
62
64
65using namespace RooStats;
66using namespace std;
67
68////////////////////////////////////////////////////////////////////////////////
69
70SPlot::~SPlot()
71{
72 if(TestBit(kOwnData) && fSData)
73 delete fSData;
74
75}
76
77////////////////////////////////////////////////////////////////////////////////
78/// Default constructor
79
81 TNamed()
82{
83 RooArgList Args;
84
85 fSWeightVars = Args;
86
87 fSData = NULL;
88
89}
90
91////////////////////////////////////////////////////////////////////////////////
92
93SPlot::SPlot(const char* name, const char* title):
94 TNamed(name, title)
95{
96 RooArgList Args;
97
98 fSWeightVars = Args;
99
100 fSData = NULL;
101
102}
103
104////////////////////////////////////////////////////////////////////////////////
105///Constructor from a RooDataSet
106///No sWeighted variables are present
107
108SPlot::SPlot(const char* name, const char* title, const RooDataSet &data):
109 TNamed(name, title)
110{
111 RooArgList Args;
112
113 fSWeightVars = Args;
114
115 fSData = (RooDataSet*) &data;
116}
117
118////////////////////////////////////////////////////////////////////////////////
119/// Copy Constructor from another SPlot
120
121SPlot::SPlot(const SPlot &other):
122 TNamed(other)
123{
124 RooArgList Args = (RooArgList) other.GetSWeightVars();
125
127
128 fSData = (RooDataSet*) other.GetSDataSet();
129
130}
131
132////////////////////////////////////////////////////////////////////////////////
133
134SPlot::SPlot(const char* name, const char* title, RooDataSet& data, RooAbsPdf* pdf,
135 const RooArgList &yieldsList, const RooArgSet &projDeps,
136 bool includeWeights, bool cloneData, const char* newName):
137 TNamed(name, title)
138{
139 if(cloneData == 1) {
140 fSData = (RooDataSet*) data.Clone(newName);
142 }
143 else
144 fSData = (RooDataSet*) &data;
145
146 // Add check that yieldsList contains all RooRealVars
147 TIterator* iter = yieldsList.createIterator() ;
148 RooAbsArg* arg ;
149 while((arg=(RooAbsArg*)iter->Next())) {
150 if (!dynamic_cast<RooRealVar*>(arg)) {
151 coutE(InputArguments) << "SPlot::SPlot(" << GetName() << ") input argument "
152 << arg->GetName() << " is not of type RooRealVar " << endl ;
153 throw string(Form("SPlot::SPlot(%s) input argument %s is not of type RooRealVar",GetName(),arg->GetName())) ;
154 }
155 }
156 delete iter ;
157
158 //Construct a new SPlot class,
159 //calculate sWeights, and include them
160 //in the RooDataSet of this class.
161
162 this->AddSWeight(pdf, yieldsList, projDeps, includeWeights);
163}
164
165////////////////////////////////////////////////////////////////////////////////
166
168{
169 if(data) {
171 return fSData;
172 } else
173 return NULL;
174}
175
176////////////////////////////////////////////////////////////////////////////////
177
179{
180 return fSData;
181}
182
183////////////////////////////////////////////////////////////////////////////////
184
185Double_t SPlot::GetSWeight(Int_t numEvent, const char* sVariable) const
186{
187 if(numEvent > fSData->numEntries() )
188 {
189 coutE(InputArguments) << "Invalid Entry Number" << endl;
190 return -1;
191 }
192
193 if(numEvent < 0)
194 {
195 coutE(InputArguments) << "Invalid Entry Number" << endl;
196 return -1;
197 }
198
199 Double_t totalYield = 0;
200
201 std::string varname(sVariable);
202 varname += "_sw";
203
204
205 if(fSWeightVars.find(sVariable) )
206 {
207 RooArgSet Row(*fSData->get(numEvent));
208 totalYield += Row.getRealValue(sVariable);
209
210 return totalYield;
211 }
212
213 if( fSWeightVars.find(varname.c_str()) )
214 {
215
216 RooArgSet Row(*fSData->get(numEvent));
217 totalYield += Row.getRealValue(varname.c_str() );
218
219 return totalYield;
220 }
221
222 else
223 coutE(InputArguments) << "InputVariable not in list of sWeighted variables" << endl;
224
225 return -1;
226}
227
228
229////////////////////////////////////////////////////////////////////////////////
230/// Sum the SWeights for a particular event.
231/// This sum should equal the total weight of that event.
232/// This method is intended to be used as a check.
233
235{
236 if(numEvent > fSData->numEntries() )
237 {
238 coutE(InputArguments) << "Invalid Entry Number" << endl;
239 return -1;
240 }
241
242 if(numEvent < 0)
243 {
244 coutE(InputArguments) << "Invalid Entry Number" << endl;
245 return -1;
246 }
247
248 Int_t numSWeightVars = this->GetNumSWeightVars();
249
250 Double_t eventSWeight = 0;
251
252 RooArgSet Row(*fSData->get(numEvent));
253
254 for (Int_t i = 0; i < numSWeightVars; i++)
255 eventSWeight += Row.getRealValue(fSWeightVars.at(i)->GetName() );
256
257 return eventSWeight;
258}
259
260////////////////////////////////////////////////////////////////////////////////
261/// Sum the SWeights for a particular specie over all events
262/// This should equal the total (weighted) yield of that specie
263/// This method is intended as a check.
264
265Double_t SPlot::GetYieldFromSWeight(const char* sVariable) const
266{
267
268 Double_t totalYield = 0;
269
270 std::string varname(sVariable);
271 varname += "_sw";
272
273
274 if(fSWeightVars.find(sVariable) )
275 {
276 for(Int_t i=0; i < fSData->numEntries(); i++)
277 {
278 RooArgSet Row(*fSData->get(i));
279 totalYield += Row.getRealValue(sVariable);
280 }
281
282 return totalYield;
283 }
284
285 if( fSWeightVars.find(varname.c_str()) )
286 {
287 for(Int_t i=0; i < fSData->numEntries(); i++)
288 {
289 RooArgSet Row(*fSData->get(i));
290 totalYield += Row.getRealValue(varname.c_str() );
291 }
292
293 return totalYield;
294 }
295
296 else
297 coutE(InputArguments) << "InputVariable not in list of sWeighted variables" << endl;
298
299 return -1;
300}
301
302
303////////////////////////////////////////////////////////////////////////////////
304/// Return a RooArgList containing the SWeights
305
307{
308
310
311 return Args;
312
313}
314
315////////////////////////////////////////////////////////////////////////////////
316/// Return the number of SWeights
317/// In other words, return the number of
318/// species that we are trying to extract.
319
321{
323
324 return Args.getSize();
325}
326
327////////////////////////////////////////////////////////////////////////////////
328/// Method which adds the sWeights to the dataset.
329/// Input is the PDF, a RooArgList of the yields (floating)
330/// and a RooArgSet of the projDeps.
331///
332/// The projDeps will not be normalized over when calculating the SWeights
333/// and will be considered parameters, not observables.
334///
335/// The SPlot will contain two new variables for each specie of name "varname":
336///
337/// L_varname is the value of the pdf for the variable "varname" at values of this event
338/// varname_sw is the value of the sWeight for the variable "varname" for this event
339///
340/// Find Parameters in the PDF to be considered fixed when calculating the SWeights
341/// and be sure to NOT include the yields in that list
342
343void SPlot::AddSWeight( RooAbsPdf* pdf, const RooArgList &yieldsTmp,
344 const RooArgSet &projDeps, bool includeWeights)
345{
346
348
349 RooArgList* constParameters = (RooArgList*)pdf->getParameters(fSData) ;
350 constParameters->remove(yieldsTmp, kTRUE, kTRUE);
351
352
353 // Set these parameters constant and store them so they can later
354 // be set to not constant
355 std::vector<RooRealVar*> constVarHolder;
356
357 for(Int_t i = 0; i < constParameters->getSize(); i++)
358 {
359 RooRealVar* varTemp = ( dynamic_cast<RooRealVar*>( constParameters->at(i) ) );
360 if(varTemp && varTemp->isConstant() == 0 )
361 {
362 varTemp->setConstant();
363 constVarHolder.push_back(varTemp);
364 }
365 }
366
367 // Fit yields to the data with all other variables held constant
368 // This is necessary because SPlot assumes the yields minimise -Log(likelihood)
369
371
372 // Hold the value of the fitted yields
373 std::vector<double> yieldsHolder;
374
375 for(Int_t i = 0; i < yieldsTmp.getSize(); i++)
376 yieldsHolder.push_back( ((RooRealVar*) yieldsTmp.at(i))->getVal());
377
378 Int_t nspec = yieldsTmp.getSize();
379 RooArgList yields = *(RooArgList*)yieldsTmp.snapshot(kFALSE);
380
381 if(currentLevel <= RooFit::DEBUG)
382 {
383 coutI(InputArguments) << "Printing Yields" << endl;
384 yields.Print();
385 }
386
387 // The list of variables to normalize over when calculating PDF values.
388
389 RooArgSet vars(*fSData->get() );
390 vars.remove(projDeps, kTRUE, kTRUE);
391
392 // Attach data set
393
394 // const_cast<RooAbsPdf*>(pdf)->attachDataSet(*fSData);
395
396 pdf->attachDataSet(*fSData);
397
398 // first calculate the pdf values for all species and all events
399 std::vector<RooRealVar*> yieldvars ;
400 RooArgSet* parameters = pdf->getParameters(fSData) ;
401
402 std::vector<Double_t> yieldvalues ;
403 for (Int_t k = 0; k < nspec; ++k)
404 {
405 RooRealVar* thisyield = dynamic_cast<RooRealVar*>(yields.at(k)) ;
406 if (thisyield) {
407 RooRealVar* yieldinpdf = dynamic_cast<RooRealVar*>(parameters->find(thisyield->GetName() )) ;
408
409 if (yieldinpdf) {
410 coutI(InputArguments)<< "yield in pdf: " << yieldinpdf->GetName() << " " << thisyield->getVal() << endl;
411
412 yieldvars.push_back(yieldinpdf) ;
413 yieldvalues.push_back(thisyield->getVal()) ;
414 }
415 }
416 }
417
418 Int_t numevents = fSData->numEntries() ;
419
420 std::vector<std::vector<Double_t> > pdfvalues(numevents,std::vector<Double_t>(nspec,0)) ;
421
422
423 // set all yield to zero
424 for(Int_t m=0; m<nspec; ++m) yieldvars[m]->setVal(0) ;
425
426
427 // For every event and for every specie,
428 // calculate the value of the component pdf for that specie
429 // by setting the yield of that specie to 1
430 // and all others to 0. Evaluate the pdf for each event
431 // and store the values.
432
433 RooArgSet * pdfvars = pdf->getVariables();
434
435 for (Int_t ievt = 0; ievt <numevents; ievt++)
436 {
437 // if (ievt % 100 == 0)
438 // coutP(Eval) << ".";
439
440
441 //FIX THIS PART, EVALUATION PROGRESS!!
442
443 RooStats::SetParameters(fSData->get(ievt), pdfvars);
444
445 // RooArgSet row(*fSData->get(ievt));
446
447 for(Int_t k = 0; k < nspec; ++k)
448 {
449 //Check that range of yields is at least (0,1), and fix otherwise
450 if(yieldvars[k]->getMin() > 0)
451 {
452 coutW(InputArguments) << "Minimum Range for " << yieldvars[k]->GetName() << " must be 0. ";
453 coutW(InputArguments) << "Setting min range to 0" << std::endl;
454 yieldvars[k]->setMin(0);
455 }
456
457 if(yieldvars[k]->getMax() < 1)
458 {
459 coutW(InputArguments) << "Maximum Range for " << yieldvars[k]->GetName() << " must be 1. ";
460 coutW(InputArguments) << "Setting max range to 1" << std::endl;
461 yieldvars[k]->setMax(1);
462 }
463
464 // set this yield to 1
465 yieldvars[k]->setVal( 1 ) ;
466 // evaluate the pdf
467 Double_t f_k = pdf->getVal(&vars) ;
468 pdfvalues[ievt][k] = f_k ;
469 if( !(f_k>1 || f_k<1) )
470 coutW(InputArguments) << "Strange pdf value: " << ievt << " " << k << " " << f_k << std::endl ;
471 yieldvars[k]->setVal( 0 ) ;
472 }
473 }
474 delete pdfvars;
475
476 // check that the likelihood normalization is fine
477 std::vector<Double_t> norm(nspec,0) ;
478 for (Int_t ievt = 0; ievt <numevents ; ievt++)
479 {
480 Double_t dnorm(0) ;
481 for(Int_t k=0; k<nspec; ++k) dnorm += yieldvalues[k] * pdfvalues[ievt][k] ;
482 for(Int_t j=0; j<nspec; ++j) norm[j] += pdfvalues[ievt][j]/dnorm ;
483 }
484
485 coutI(Contents) << "likelihood norms: " ;
486
487 for(Int_t k=0; k<nspec; ++k) coutI(Contents) << norm[k] << " " ;
488 coutI(Contents) << std::endl ;
489
490 // Make a TMatrixD to hold the covariance matrix.
491 TMatrixD covInv(nspec, nspec);
492 for (Int_t i = 0; i < nspec; i++) for (Int_t j = 0; j < nspec; j++) covInv(i,j) = 0;
493
494 coutI(Contents) << "Calculating covariance matrix";
495
496
497 // Calculate the inverse covariance matrix, using weights
498 for (Int_t ievt = 0; ievt < numevents; ++ievt)
499 {
500
501 fSData->get(ievt) ;
502
503 // Calculate contribution to the inverse of the covariance
504 // matrix. See BAD 509 V2 eqn. 15
505
506 // Sum for the denominator
507 Double_t dsum(0);
508 for(Int_t k = 0; k < nspec; ++k)
509 dsum += pdfvalues[ievt][k] * yieldvalues[k] ;
510
511 for(Int_t n=0; n<nspec; ++n)
512 for(Int_t j=0; j<nspec; ++j)
513 {
514 if(includeWeights == kTRUE)
515 covInv(n,j) += fSData->weight()*pdfvalues[ievt][n]*pdfvalues[ievt][j]/(dsum*dsum) ;
516 else
517 covInv(n,j) += pdfvalues[ievt][n]*pdfvalues[ievt][j]/(dsum*dsum) ;
518 }
519
520 //ADDED WEIGHT ABOVE
521
522 }
523
524 // Covariance inverse should now be computed!
525
526 // Invert to get the covariance matrix
527 if (covInv.Determinant() <=0)
528 {
529 coutE(Eval) << "SPlot Error: covariance matrix is singular; I can't invert it!" << std::endl;
530 covInv.Print();
531 return;
532 }
533
534 TMatrixD covMatrix(TMatrixD::kInverted,covInv);
535
536 //check cov normalization
537 if(currentLevel <= RooFit::DEBUG)
538 {
539 coutI(Eval) << "Checking Likelihood normalization: " << std::endl;
540 coutI(Eval) << "Yield of specie Sum of Row in Matrix Norm" << std::endl;
541 for(Int_t k=0; k<nspec; ++k)
542 {
543 Double_t covnorm(0) ;
544 for(Int_t m=0; m<nspec; ++m) covnorm += covInv[k][m]*yieldvalues[m] ;
545 Double_t sumrow(0) ;
546 for(Int_t m = 0; m < nspec; ++m) sumrow += covMatrix[k][m] ;
547 coutI(Eval) << yieldvalues[k] << " " << sumrow << " " << covnorm << endl ;
548 }
549 }
550
551 // calculate for each event the sWeight (BAD 509 V2 eq. 21)
552 coutI(Eval) << "Calculating sWeight" << std::endl;
553 std::vector<RooRealVar*> sweightvec ;
554 std::vector<RooRealVar*> pdfvec ;
555 RooArgSet sweightset ;
556
557 // Create and label the variables
558 // used to store the SWeights
559
561
562 for(Int_t k=0; k<nspec; ++k)
563 {
564 std::string wname = std::string(yieldvars[k]->GetName()) + "_sw";
565 RooRealVar* var = new RooRealVar(wname.c_str(),wname.c_str(),0) ;
566 sweightvec.push_back( var) ;
567 sweightset.add(*var) ;
568 fSWeightVars.add(*var);
569
570 wname = "L_" + std::string(yieldvars[k]->GetName());
571 var = new RooRealVar(wname.c_str(),wname.c_str(),0) ;
572 pdfvec.push_back( var) ;
573 sweightset.add(*var) ;
574 }
575
576 // Create and fill a RooDataSet
577 // with the SWeights
578
579 RooDataSet* sWeightData = new RooDataSet("dataset", "dataset with sWeights", sweightset);
580
581 for(Int_t ievt = 0; ievt < numevents; ++ievt)
582 {
583
584 fSData->get(ievt) ;
585
586 // sum for denominator
587 Double_t dsum(0);
588 for(Int_t k = 0; k < nspec; ++k) dsum += pdfvalues[ievt][k] * yieldvalues[k] ;
589 // covariance weighted pdf for each specief
590 for(Int_t n=0; n<nspec; ++n)
591 {
592 Double_t nsum(0) ;
593 for(Int_t j=0; j<nspec; ++j) nsum += covMatrix(n,j) * pdfvalues[ievt][j] ;
594
595
596 //Add the sWeights here!!
597 //Include weights,
598 //ie events weights are absorbed into sWeight
599
600
601 if(includeWeights == kTRUE) sweightvec[n]->setVal(fSData->weight() * nsum/dsum) ;
602 else sweightvec[n]->setVal( nsum/dsum) ;
603
604 pdfvec[n]->setVal( pdfvalues[ievt][n] ) ;
605
606 if( !(fabs(nsum/dsum)>=0 ) )
607 {
608 coutE(Contents) << "error: " << nsum/dsum << endl ;
609 return;
610 }
611 }
612
613 sWeightData->add(sweightset) ;
614 }
615
616
617 // Add the SWeights to the original data set
618
619 fSData->merge(sWeightData);
620
621 delete sWeightData;
622
623 //Restore yield values
624
625 for(Int_t i = 0; i < yieldsTmp.getSize(); i++)
626 ((RooRealVar*) yieldsTmp.at(i))->setVal(yieldsHolder.at(i));
627
628 //Make any variables that were forced to constant no longer constant
629
630 for(Int_t i=0; i < (Int_t) constVarHolder.size(); i++)
631 constVarHolder.at(i)->setConstant(kFALSE);
632
633 return;
634
635}
#define coutI(a)
Definition: RooMsgService.h:31
#define coutW(a)
Definition: RooMsgService.h:33
#define coutE(a)
Definition: RooMsgService.h:34
int Int_t
Definition: RtypesCore.h:41
const Bool_t kFALSE
Definition: RtypesCore.h:88
double Double_t
Definition: RtypesCore.h:55
const Bool_t kTRUE
Definition: RtypesCore.h:87
#define ClassImp(name)
Definition: Rtypes.h:363
char * Form(const char *fmt,...)
RooAbsArg is the common abstract base class for objects that represent a value (of arbitrary type) an...
Definition: RooAbsArg.h:66
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:533
RooArgSet * getVariables(Bool_t stripDisconnected=kTRUE) const
Return RooArgSet with all variables (tree leaf nodes of expresssion tree)
Definition: RooAbsArg.cxx:2074
Bool_t isConstant() const
Definition: RooAbsArg.h:266
void attachDataSet(const RooAbsData &set)
Replace server nodes with names matching the dataset variable names with those data set variables,...
Definition: RooAbsArg.cxx:1489
Int_t getSize() const
RooAbsCollection * snapshot(Bool_t deepCopy=kTRUE) const
Take a snap shot of current collection contents: An owning collection is returned containing clones o...
virtual RooAbsArg * addClone(const RooAbsArg &var, Bool_t silent=kFALSE)
Add a clone of the specified argument to list.
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
virtual void Print(Option_t *options=0) const
This method must be overridden when a class wants to print itself.
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE)
Remove the specified argument from our list.
TIterator * createIterator(Bool_t dir=kIterForward) const
RooAbsArg * find(const char *name) const
Find object with given name in list.
virtual Int_t numEntries() const
Definition: RooAbsData.cxx:285
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Definition: RooAbsPdf.h:41
virtual RooFitResult * fitTo(RooAbsData &data, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none())
Fit PDF to given dataset.
Definition: RooAbsPdf.cxx:1081
void setConstant(Bool_t value=kTRUE)
Double_t getVal(const RooArgSet *set=0) const
Evaluate object. Returns either cached value or triggers a recalculation.
Definition: RooAbsReal.h:64
RooAbsArg * at(Int_t idx) const
Definition: RooArgList.h:84
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:28
Double_t getRealValue(const char *name, Double_t defVal=0, Bool_t verbose=kFALSE) const
Get value of a RooAbsReal stored in set with given name.
Definition: RooArgSet.cxx:472
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
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:31
virtual const RooArgSet * get(Int_t index) const
Return RooArgSet with coordinates of event 'index'.
Definition: RooDataSet.cxx:995
virtual Double_t weight() const
Return event weight of current event.
Definition: RooDataSet.cxx:955
virtual void add(const RooArgSet &row, Double_t weight=1.0, Double_t weightError=0)
Add a data point, with its coordinates specified in the 'data' argset, to the data set.
Bool_t merge(RooDataSet *data1, RooDataSet *data2=0, RooDataSet *data3=0, RooDataSet *data4=0, RooDataSet *data5=0, RooDataSet *data6=0)
static RooMsgService & instance()
Return reference to singleton instance.
RooFit::MsgLevel globalKillBelow() const
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:36
This class calculates sWeights used to create an sPlot.
Definition: SPlot.h:32
Double_t GetYieldFromSWeight(const char *sVariable) const
Sum the SWeights for a particular specie over all events This should equal the total (weighted) yield...
Definition: SPlot.cxx:265
Double_t GetSWeight(Int_t numEvent, const char *sVariable) const
Definition: SPlot.cxx:185
RooDataSet * fSData
Definition: SPlot.h:74
RooArgList GetSWeightVars() const
Return a RooArgList containing the SWeights.
Definition: SPlot.cxx:306
SPlot()
Default constructor.
Definition: SPlot.cxx:80
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:320
void AddSWeight(RooAbsPdf *pdf, const RooArgList &yieldsTmp, const RooArgSet &projDeps=RooArgSet(), bool includeWeights=kTRUE)
Method which adds the sWeights to the dataset.
Definition: SPlot.cxx:343
RooDataSet * SetSData(RooDataSet *data)
Definition: SPlot.cxx:167
RooDataSet * GetSDataSet() const
Definition: SPlot.cxx:178
RooArgList fSWeightVars
Definition: SPlot.h:70
Double_t GetSumOfEventSWeight(Int_t numEvent) const
Sum the SWeights for a particular event.
Definition: SPlot.cxx:234
Iterator abstract base class.
Definition: TIterator.h:30
virtual TObject * Next()=0
void Print(Option_t *name="") const
Print the matrix as a table of elements.
virtual Double_t Determinant() const
Return the matrix determinant.
Definition: TMatrixT.cxx:1361
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
virtual void Clear(Option_t *="")
Definition: TObject.h:100
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:172
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:694
const Int_t n
Definition: legend1.C:16
VecExpr< UnaryOp< Fabs< T >, VecExpr< A, T, D >, T >, T, D > fabs(const VecExpr< A, T, D > &rhs)
RooCmdArg SumW2Error(Bool_t flag)
RooCmdArg Extended(Bool_t flag=kTRUE)
RooCmdArg PrintEvalErrors(Int_t numErrors)
@ InputArguments
Definition: RooGlobalFunc.h:58
RooCmdArg PrintLevel(Int_t code)
@(#)root/roostats:$Id$ Author: George Lewis, Kyle Cranmer
Definition: Asimov.h:20
void SetParameters(const RooArgSet *desiredVals, RooArgSet *paramsToChange)
Definition: RooStatsUtils.h:58
STL namespace.
auto * m
Definition: textangle.C:8