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 <iostream>
19#include <cmath>
20
21#include "RooAbsReal.h"
22#include "RooRealVar.h"
23#include "RooArgList.h"
24#include "RooMsgService.h"
25
27
28using std::vector, std::endl;
29
31
32using namespace RooStats;
33using namespace HistFactory;
34
35////////////////////////////////////////////////////////////////////////////////
36
37LinInterpVar::LinInterpVar(const char* name, const char* title,
38 const RooArgList& paramList,
39 double nominal, vector<double> low, vector<double> high) :
40 RooAbsReal(name, title),
41 _paramList("paramList","List of paramficients",this),
42 _nominal(nominal), _low(low), _high(high)
43{
44
45 for (auto param : paramList) {
46 if (!dynamic_cast<RooAbsReal*>(param)) {
47 coutE(InputArguments) << "LinInterpVar::ctor(" << GetName() << ") ERROR: paramficient " << param->GetName()
48 << " is not of type RooAbsReal" << endl ;
49 assert(0) ;
50 }
51 _paramList.add(*param) ;
52 }
53}
54
55////////////////////////////////////////////////////////////////////////////////
56/// Constructor of flat polynomial function
57
58LinInterpVar::LinInterpVar(const char* name, const char* title) :
59 RooAbsReal(name, title),
60 _paramList("paramList","List of coefficients",this),
61 _nominal(0)
62{
63
64}
65
66////////////////////////////////////////////////////////////////////////////////
67
68LinInterpVar::LinInterpVar(const LinInterpVar& other, const char* name) :
69 RooAbsReal(other, name),
70 _paramList("paramList",this,other._paramList),
71 _nominal(other._nominal), _low(other._low), _high(other._high)
72
73{
74 // Copy constructor
75
76}
77
78
79////////////////////////////////////////////////////////////////////////////////
80/// Calculate and return value of polynomial
81
83{
84 double sum(_nominal) ;
85
86 int i=0;
87 for(auto const* param: static_range_cast<RooAbsReal *>(_paramList)) {
88 if (param->getVal() > 0) {
89 sum += param->getVal()*(_high.at(i) - _nominal );
90 } else {
91 sum += param->getVal() * (_nominal - _low.at(i));
92 }
93
94 ++i;
95 }
96
97 if(sum<=0) {
98 sum=1E-9;
99 }
100
101 return sum;
102}
103
104
105
#define coutE(a)
#define ClassImp(name)
Definition Rtypes.h:377
char name[80]
Definition TGX11.cxx:110
Abstract base class for objects that represent a real value and implements functionality common to al...
Definition RooAbsReal.h:59
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition RooArgList.h:22
bool add(const RooAbsArg &var, bool valueServer, bool shapeServer, bool silent)
Overloaded RooCollection_t::add() method insert object into set and registers object as server to own...
RooAbsReal that does piecewise-linear interpolations.
double evaluate() const override
Calculate and return value of polynomial.
const char * GetName() const override
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