Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
LinInterpVar.cxx
Go to the documentation of this file.
1// @(#)root/roostats:$Id: cranmer $
2// Author: Kyle Cranmer, Akira Shibata
3/*************************************************************************
4 * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers. *
5 * All rights reserved. *
6 * *
7 * For the licensing terms see $ROOTSYS/LICENSE. *
8 * For the list of contributors see $ROOTSYS/README/CREDITS. *
9 *************************************************************************/
10
11////////////////////////////////////////////////////////////////////////////////
12
13/** \class RooStats::HistFactory::LinInterpVar
14 * \ingroup HistFactory
15 * RooAbsReal that does piecewise-linear interpolations.
16 */
17
18#include "RooFit.h"
19
20#include <iostream>
21#include <cmath>
22
23#include "RooAbsReal.h"
24#include "RooRealVar.h"
25#include "RooArgList.h"
26#include "RooMsgService.h"
27
29
30using namespace std;
31
33
34using namespace RooStats;
35using namespace HistFactory;
36
37////////////////////////////////////////////////////////////////////////////////
38/// Default constructor
39
41{
43 _nominal = 0 ;
44}
45
46
47////////////////////////////////////////////////////////////////////////////////
48
49LinInterpVar::LinInterpVar(const char* name, const char* title,
50 const RooArgList& paramList,
51 double nominal, vector<double> low, vector<double> high) :
52 RooAbsReal(name, title),
53 _paramList("paramList","List of paramficients",this),
54 _nominal(nominal), _low(low), _high(high)
55{
57
58
59 TIterator* paramIter = paramList.createIterator() ;
60 RooAbsArg* param ;
61 while((param = (RooAbsArg*)paramIter->Next())) {
62 if (!dynamic_cast<RooAbsReal*>(param)) {
63 coutE(InputArguments) << "LinInterpVar::ctor(" << GetName() << ") ERROR: paramficient " << param->GetName()
64 << " is not of type RooAbsReal" << endl ;
65 assert(0) ;
66 }
67 _paramList.add(*param) ;
68 }
69 delete paramIter ;
70
71}
72
73////////////////////////////////////////////////////////////////////////////////
74/// Constructor of flat polynomial function
75
76LinInterpVar::LinInterpVar(const char* name, const char* title) :
77 RooAbsReal(name, title),
78 _paramList("paramList","List of coefficients",this),
79 _nominal(0)
80{
82}
83
84////////////////////////////////////////////////////////////////////////////////
85
86LinInterpVar::LinInterpVar(const LinInterpVar& other, const char* name) :
87 RooAbsReal(other, name),
88 _paramList("paramList",this,other._paramList),
89 _nominal(other._nominal), _low(other._low), _high(other._high)
90
91{
92 // Copy constructor
94
95}
96
97
98////////////////////////////////////////////////////////////////////////////////
99/// Destructor
100
102{
103 delete _paramIter ;
104}
105
106
107
108
109////////////////////////////////////////////////////////////////////////////////
110/// Calculate and return value of polynomial
111
113{
115 _paramIter->Reset() ;
116
117 RooAbsReal* param ;
118 //const RooArgSet* nset = _paramList.nset() ;
119 int i=0;
120
121 while((param=(RooAbsReal*)_paramIter->Next())) {
122 // param->Print("v");
123
124 if(param->getVal()>0)
125 sum += param->getVal()*(_high.at(i) - _nominal );
126 else
127 sum += param->getVal()*(_nominal - _low.at(i));
128
129 ++i;
130 }
131
132 if(sum<=0) {
133 sum=1E-9;
134 }
135
136 return sum;
137}
138
139
140
#define coutE(a)
#define ClassImp(name)
Definition Rtypes.h:364
char name[80]
Definition TGX11.cxx:110
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition RooAbsArg.h:72
TIterator * createIterator(Bool_t dir=kIterForward) const
TIterator-style iteration over contained elements.
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition RooAbsReal.h:61
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition RooAbsReal.h:91
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition RooArgList.h:21
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE) override
Reimplementation of standard RooArgList::add()
RooAbsReal that does piecewise-linear interpolations.
Double_t evaluate() const
do not persist
Iterator abstract base class.
Definition TIterator.h:30
virtual void Reset()=0
virtual TObject * Next()=0
virtual const char * GetName() const
Returns name of object.
Definition TNamed.h:47
Namespace for the RooStats classes.
Definition Asimov.h:19
static uint64_t sum(uint64_t i)
Definition Factory.cxx:2345