MINT2
Public Member Functions | Private Member Functions | Private Attributes | Friends | List of all members
RooCubicSplineFun Class Reference

#include <RooCubicSplineFun.h>

Inheritance diagram for RooCubicSplineFun:
RooAbsGaussModelEfficiency

Public Member Functions

 RooCubicSplineFun ()
 
 RooCubicSplineFun (const char *name, const char *title, RooRealVar &x, const std::vector< double > &knots, const std::vector< double > &values, const std::vector< double > &errors=std::vector< double >(), double smooth=0, bool constCoeffs=true)
 
 RooCubicSplineFun (const char *name, const char *title, RooRealVar &x, const TGraph *graph, bool constCoeffs=true)
 
 RooCubicSplineFun (const char *name, const char *title, RooRealVar &x, const TH1 *hist, double smooth=0, bool constCoeffs=true)
 
 RooCubicSplineFun (const char *name, const char *title, RooRealVar &x, const TGraphErrors *graph, double smooth=0, bool constCoeffs=true)
 
 RooCubicSplineFun (const char *name, const char *title, RooRealVar &x, const char *knotBinningName, const RooArgList &coefList)
 
 RooCubicSplineFun (const char *name, const char *title, RooRealVar &x, const std::vector< double > &knots, const RooArgList &coefList)
 
 ~RooCubicSplineFun ()
 
 RooCubicSplineFun (const RooCubicSplineFun &other, const char *name=0)
 
TObject * clone (const char *newname) const
 
Int_t getAnalyticalIntegral (RooArgSet &allVars, RooArgSet &analVars, const char *rangeName) const
 
Double_t analyticalIntegral (Int_t code, const char *rangeName) const
 
Int_t getMaxVal (const RooArgSet &vars) const
 
Double_t maxVal (Int_t code) const
 
std::complex< double > productAnalyticalIntegral (Double_t umin, Double_t umax, Double_t scale, Double_t offset, const std::complex< double > &z) const
 
unsigned knotSize () const
 
double u (int i) const
 
const std::vector< double > & knots () const
 
const RooArgList & coefficients () const
 
- Public Member Functions inherited from RooAbsGaussModelEfficiency
 RooAbsGaussModelEfficiency ()
 
 RooAbsGaussModelEfficiency (const char *name, const char *title, const char *unit="")
 
 RooAbsGaussModelEfficiency (const RooAbsGaussModelEfficiency &other, const char *name=0)
 
 ~RooAbsGaussModelEfficiency ()
 

Private Member Functions

void init (const char *name, const std::vector< double > &heights, const std::vector< double > &errors, double smooth, bool constCoeffs)
 
Double_t evaluate () const
 
std::complex< double > gaussIntegralE (bool left, const RooGaussModelAcceptance::M_n< 4U > &dM, const RooGaussModelAcceptance::K_n &K, double offset, double *sc) const
 

Private Attributes

RooRealProxy _x
 
RooListProxy _coefList
 
RooCubicSplineKnot _aux
 

Friends

class RooSplineProduct
 

Detailed Description

Definition at line 27 of file RooCubicSplineFun.h.

Constructor & Destructor Documentation

◆ RooCubicSplineFun() [1/8]

RooCubicSplineFun::RooCubicSplineFun ( )

Definition at line 79 of file RooCubicSplineFun.cpp.

80  : _aux(0,0)
81 {
82 }
RooCubicSplineKnot _aux

◆ RooCubicSplineFun() [2/8]

RooCubicSplineFun::RooCubicSplineFun ( const char *  name,
const char *  title,
RooRealVar &  x,
const std::vector< double > &  knots,
const std::vector< double > &  values,
const std::vector< double > &  errors = std::vector< double >(),
double  smooth = 0,
bool  constCoeffs = true 
)

◆ RooCubicSplineFun() [3/8]

RooCubicSplineFun::RooCubicSplineFun ( const char *  name,
const char *  title,
RooRealVar &  x,
const TGraph *  graph,
bool  constCoeffs = true 
)

Definition at line 99 of file RooCubicSplineFun.cpp.

100  :
101  RooAbsGaussModelEfficiency(name, title),
102  _x("x", "Dependent", this, x),
103  _coefList("coefficients","List of coefficients",this),
104  _aux(0,0)
105 {
106  int nPoints = graph->GetN();
107  vector<double> centres, values;
108  for (int i=0;i<nPoints ;++i) {
109  double x,y;
110  graph->GetPoint(i,x,y);
111  centres.push_back(x);
112  values.push_back(y);
113  }
114  _aux = RooCubicSplineKnot( centres.begin(), centres.end() );
115  vector<double> errs;
116  init(name, values, errs, 0, constCoeffs);
117 }
void init(const char *name, const std::vector< double > &heights, const std::vector< double > &errors, double smooth, bool constCoeffs)
RooCubicSplineKnot _aux
RooListProxy _coefList

◆ RooCubicSplineFun() [4/8]

RooCubicSplineFun::RooCubicSplineFun ( const char *  name,
const char *  title,
RooRealVar &  x,
const TH1 *  hist,
double  smooth = 0,
bool  constCoeffs = true 
)

Definition at line 120 of file RooCubicSplineFun.cpp.

122  :
123  RooAbsGaussModelEfficiency(name, title),
124  _x("x", "Dependent", this, x),
125  _coefList("coefficients","List of coefficients",this),
126  _aux(0,0)
127 {
128  // bin 0 is underflow, and bin nBins + 1 is overflow...
129  int nBins = hist->GetNbinsX();
130  vector<double> centres, values;
131  for (int i=0;i<nBins ;++i) {
132  centres.push_back(hist->GetBinCenter(1+i));
133  values.push_back(hist->GetBinContent(1+i));
134  }
135  _aux = RooCubicSplineKnot( centres.begin(), centres.end() );
136 
137  vector<double> errs;
138  if (smooth>0) for (int i=0;i<nBins ;++i) errs.push_back(hist->GetBinError(1+i));
139 
140  init(name, values, errs, smooth, constCoeffs);
141 }
void init(const char *name, const std::vector< double > &heights, const std::vector< double > &errors, double smooth, bool constCoeffs)
RooCubicSplineKnot _aux
RooListProxy _coefList

◆ RooCubicSplineFun() [5/8]

RooCubicSplineFun::RooCubicSplineFun ( const char *  name,
const char *  title,
RooRealVar &  x,
const TGraphErrors *  graph,
double  smooth = 0,
bool  constCoeffs = true 
)

Definition at line 143 of file RooCubicSplineFun.cpp.

145  :
146  RooAbsGaussModelEfficiency(name, title),
147  _x("x", "Dependent", this, x),
148  _coefList("coefficients","List of coefficients",this),
149  _aux(0,0)
150 {
151  int nPoints = graph->GetN();
152  vector<double> centres, values;
153  for (int i=0;i<nPoints ;++i) {
154  double x,y;
155  graph->GetPoint(i,x,y);
156  centres.push_back(x);
157  values.push_back(y);
158  }
159  _aux = RooCubicSplineKnot( centres.begin(), centres.end() );
160 
161  vector<double> errs;
162  if (smooth>0) for (int i=0;i<nPoints ;++i) errs.push_back(graph->GetErrorY(i));
163 
164  init(name, values, errs, smooth, constCoeffs);
165 }
void init(const char *name, const std::vector< double > &heights, const std::vector< double > &errors, double smooth, bool constCoeffs)
RooCubicSplineKnot _aux
RooListProxy _coefList

◆ RooCubicSplineFun() [6/8]

RooCubicSplineFun::RooCubicSplineFun ( const char *  name,
const char *  title,
RooRealVar &  x,
const char *  knotBinningName,
const RooArgList &  coefList 
)

Definition at line 168 of file RooCubicSplineFun.cpp.

170  :
171  RooAbsGaussModelEfficiency(name, title),
172  _x("x", "Dependent", this, x),
173  _coefList("coefficients", "List of coefficients", this),
174  _aux(0,0)
175 {
176  // TODO: verify coefList is consistent with knots as specified by the knotBinningName binning
177  // should be N+2 coefficients for N bins...
178  const RooAbsBinning* binning = x.getBinningPtr(knotBinningName);
179  assert( binning!=0);
180  if ( coefList.getSize()!=2+binning->numBoundaries()) {
181  cout << TString::Format( "you have specified %d coefficients for %d knots. The differentce should 2!",coefList.getSize(),binning->numBoundaries()) << endl;
182  throw TString::Format( "you have specified %d coefficients for %d knots. The differentce should 2!",coefList.getSize(),binning->numBoundaries());
183  }
184  _coefList.add(coefList);
185 
186  Double_t* boundaries = binning->array();
187  _aux = RooCubicSplineKnot( boundaries, boundaries + binning->numBoundaries() );
188 }
RooCubicSplineKnot _aux
RooListProxy _coefList

◆ RooCubicSplineFun() [7/8]

RooCubicSplineFun::RooCubicSplineFun ( const char *  name,
const char *  title,
RooRealVar &  x,
const std::vector< double > &  knots,
const RooArgList &  coefList 
)

◆ ~RooCubicSplineFun()

RooCubicSplineFun::~RooCubicSplineFun ( )

Definition at line 213 of file RooCubicSplineFun.cpp.

214 {
215 }

◆ RooCubicSplineFun() [8/8]

RooCubicSplineFun::RooCubicSplineFun ( const RooCubicSplineFun other,
const char *  name = 0 
)

Definition at line 204 of file RooCubicSplineFun.cpp.

204  :
205  RooAbsGaussModelEfficiency(other, name),
206  _x("x", this, other._x),
207  _coefList("coefList", this, other._coefList),
208  _aux(other._aux)
209 {
210 }
RooCubicSplineKnot _aux
RooListProxy _coefList

Member Function Documentation

◆ analyticalIntegral()

Double_t RooCubicSplineFun::analyticalIntegral ( Int_t  code,
const char *  rangeName 
) const

Definition at line 232 of file RooCubicSplineFun.cpp.

233 {
234  if (code != 1) {
235  coutE(InputArguments) << "RooCubicSplineFun::analyticalIntegral(" << GetName()
236  << "): argument \"code\" can only have value 1" << std::endl;
237  assert(code==1) ;
238  }
240 }
double analyticalIntegral(const RooArgList &b) const
RooCubicSplineKnot _aux
RooListProxy _coefList

◆ clone()

TObject* RooCubicSplineFun::clone ( const char *  newname) const
inline

Definition at line 53 of file RooCubicSplineFun.h.

53 { return new RooCubicSplineFun(*this, newname); }

◆ coefficients()

const RooArgList& RooCubicSplineFun::coefficients ( ) const
inline

Definition at line 69 of file RooCubicSplineFun.h.

69 { return _coefList; }
RooListProxy _coefList

◆ evaluate()

Double_t RooCubicSplineFun::evaluate ( ) const
private

Definition at line 218 of file RooCubicSplineFun.cpp.

219 {
220  return _aux.evaluate(_x,_coefList);
221 }
double evaluate(double _u, const RooArgList &b) const
RooCubicSplineKnot _aux
RooListProxy _coefList

◆ gaussIntegralE()

complex< double > RooCubicSplineFun::gaussIntegralE ( bool  left,
const RooGaussModelAcceptance::M_n< 4U > &  dM,
const RooGaussModelAcceptance::K_n K,
double  offset,
double *  sc 
) const
private

Definition at line 243 of file RooCubicSplineFun.cpp.

246 {
247  RooCubicSplineKnot::S_edge s_jk( _aux.S_jk_edge( left, _coefList ), offset );
248  return dM(0)*( s_jk(0,0) * K(0) * sc[0] + s_jk(0,1) * K(1) * sc[1] )
249  + dM(1)*( s_jk(1,0) * K(0) * sc[1] );
250 }
RooCubicSplineKnot::S_edge S_jk_edge(bool left, const RooArgList &b) const
RooCubicSplineKnot _aux
RooListProxy _coefList

◆ getAnalyticalIntegral()

Int_t RooCubicSplineFun::getAnalyticalIntegral ( RooArgSet &  allVars,
RooArgSet &  analVars,
const char *  rangeName 
) const

Definition at line 224 of file RooCubicSplineFun.cpp.

225 {
226  // No analytical calculation available (yet) of integrals over subranges
227  if (_x.min(rangeName)!=_aux.knots().front() || _x.max(rangeName)!=_aux.knots().back() ) return 0;
228  if (matchArgs(allVars, analVars, _x)) return 1;
229  return 0;
230 }
const std::vector< double > & knots() const
RooCubicSplineKnot _aux

◆ getMaxVal()

Int_t RooCubicSplineFun::getMaxVal ( const RooArgSet &  vars) const

Definition at line 284 of file RooCubicSplineFun.cpp.

285 {
286  // check that vars only contains _x...
287  return ( vars.getSize() == 1 && vars.contains( _x.arg() ) ) ? 1 : 0;
288 }

◆ init()

void RooCubicSplineFun::init ( const char *  name,
const std::vector< double > &  heights,
const std::vector< double > &  errors,
double  smooth,
bool  constCoeffs 
)
private

Definition at line 49 of file RooCubicSplineFun.cpp.

52  {
53  vector<double> values(heights);
54  if ( smooth > 0 ) {
55  assert(int(errors.size()) == _aux.size());
56  _aux.smooth( values, errors, smooth );
57  }
58  _aux.computeCoefficients( values );
59  for (unsigned int i=0;i<values.size();++i) {
60  if (constCoeffs) {
61  _coefList.add( RooFit::RooConst( values[i] ) );
62  } else {
63  stringstream name_str;
64  name_str << name << "_bin_" << i;
65  string n(name_str.str());
66  RooRealVar* coeff = new RooRealVar(n.c_str(), n.c_str(), values[i], 0.0, 10.0);
67 // if (i == 0 || i == values.size() - 1) coeff->setConstant(true);
68  if (i == values.size() - 2){
69  coeff->setVal(1.0);
70  coeff->setConstant(true);
71  }
72  _coefList.add(*coeff);
73  addOwnedComponents( *coeff );
74  }
75  }
76 }
void smooth(std::vector< double > &y, const std::vector< double > &dy, double lambda) const
void computeCoefficients(std::vector< double > &y, BoundaryConditions bc=BoundaryConditions()) const
RooCubicSplineKnot _aux
RooListProxy _coefList

◆ knots()

const std::vector<double>& RooCubicSplineFun::knots ( ) const
inline

Definition at line 68 of file RooCubicSplineFun.h.

68 { return _aux.knots(); }
const std::vector< double > & knots() const
RooCubicSplineKnot _aux

◆ knotSize()

unsigned RooCubicSplineFun::knotSize ( ) const
inline

Definition at line 66 of file RooCubicSplineFun.h.

66 { return _aux.size(); }
RooCubicSplineKnot _aux

◆ maxVal()

Double_t RooCubicSplineFun::maxVal ( Int_t  code) const

Definition at line 290 of file RooCubicSplineFun.cpp.

291 {
292  if (code != 1) {
293  coutE(InputArguments) << "RooCubicSplineFun::maxVal(" << GetName()
294  << "): argument \"code\" can only have value 1" << std::endl;
295  assert(code==1) ;
296  }
297  RooFIter iter = _coefList.fwdIterator();
298  RooAbsReal *c(0);
299  double res = 0;
300  while((c=(RooAbsReal*)iter.next())) {
301  double x = fabs(c->getVal());
302  if (x>res) { res = x; }
303  }
304  return res;
305 }
RooListProxy _coefList

◆ productAnalyticalIntegral()

std::complex< double > RooCubicSplineFun::productAnalyticalIntegral ( Double_t  umin,
Double_t  umax,
Double_t  scale,
Double_t  offset,
const std::complex< double > &  z 
) const
virtual

Implements RooAbsGaussModelEfficiency.

Definition at line 254 of file RooCubicSplineFun.cpp.

257 {
259  assert(knotSize()>1);
260  double lo = scale*umin+offset;
261  double hi = scale*umax+offset+1.e-7;
263  std::vector<M_n> M; M.reserve( knotSize() );
264  for (unsigned int i=0;i<knotSize();++i) {
265  double x = (u(i)-offset)/scale ;
266  // TODO: remove unnecessary calculation of integral (that gives 0) if lo>u(i) and/or u(i+1)<hi
267  if (lo>=u(i)) x = umin ; // take M[i] if lo<=u(i) else M_n(lo)
268  if (u(i)>=hi) x = umax ; // take M[i+1] if u(i+1)<=hi else M_n(hi)
269  M.push_back( M_n( x, z ) );
270  }
271  double sc[4]; for (int i=0;i<4;++i) sc[i] = pow(scale,i);
272  std::complex<double> sum(0,0);
273  if (lo<u(0)) sum += gaussIntegralE(true, M.front()-M_n( umin,z), K, offset, sc);
274  for (unsigned i=0;i<knotSize()-1 && u(i)<hi ;++i) {
275  M_n dM = M[i+1]-M[i];
276  RooCubicSplineKnot::S_jk s_jk( _aux.S_jk_sum( i, _coefList ), offset ); // pass sc into S_jk, remove from loop
277  for (int j=0;j<4;++j) for (int k=0;k<4-j;++k) sum += dM(j)*s_jk(j,k)*K(k)*sc[j+k];
278  }
279  if (hi>u(knotSize()-1)) sum += gaussIntegralE(false, M_n(umax,z)-M.back(), K, offset, sc);
280  return sum;
281 }
double u(int i) const
RooCubicSplineKnot::S_jk S_jk_sum(int i, const RooArgList &b) const
unsigned knotSize() const
std::complex< double > gaussIntegralE(bool left, const RooGaussModelAcceptance::M_n< 4U > &dM, const RooGaussModelAcceptance::K_n &K, double offset, double *sc) const
RooCubicSplineKnot _aux
RooListProxy _coefList

◆ u()

double RooCubicSplineFun::u ( int  i) const
inline

Definition at line 67 of file RooCubicSplineFun.h.

67 { return _aux.u(i); }
double u(int i) const
RooCubicSplineKnot _aux

Friends And Related Function Documentation

◆ RooSplineProduct

friend class RooSplineProduct
friend

Definition at line 28 of file RooCubicSplineFun.h.

Member Data Documentation

◆ _aux

RooCubicSplineKnot RooCubicSplineFun::_aux
private

Definition at line 75 of file RooCubicSplineFun.h.

◆ _coefList

RooListProxy RooCubicSplineFun::_coefList
private

Definition at line 74 of file RooCubicSplineFun.h.

◆ _x

RooRealProxy RooCubicSplineFun::_x
private

Definition at line 73 of file RooCubicSplineFun.h.


The documentation for this class was generated from the following files: