#include <MELANCSplineCore.h>
|
| MELANCSplineCore () |
|
| MELANCSplineCore (const char *name, const char *title) |
|
| MELANCSplineCore (const char *name, const char *title, RooAbsReal &inXVar, const std::vector< T > &inXList, Bool_t inUseFloor=true, T inFloorEval=1e-15, T inFloorInt=1e-10) |
|
| MELANCSplineCore (const MELANCSplineCore &other, const char *name=0) |
|
virtual TObject * | clone (const char *newname) const =0 |
|
virtual | ~MELANCSplineCore () |
|
virtual void | setVerbosity (VerbosityLevel flag) |
|
void | setEvalFloor (T val) |
|
void | setIntFloor (T val) |
|
void | doFloor (Bool_t flag) |
|
virtual void | setRangeValidity (const T valmin, const T valmax, const Int_t whichDirection)=0 |
|
virtual Int_t | getAnalyticalIntegral (RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const =0 |
|
virtual Double_t | analyticalIntegral (Int_t code, const char *rangeName=0) const =0 |
|
|
virtual void | emptyFcnList ()=0 |
|
void | getLeafDependents (RooRealProxy &proxy, RooArgSet &set) |
|
void | addLeafDependents (RooArgSet &set) |
|
unsigned int | npointsX () const |
|
virtual Int_t | getWhichBin (const T &val, const Int_t whichDirection) const =0 |
|
virtual void | getKappas (std::vector< T > &kappas, const Int_t whichDirection)=0 |
|
virtual T | getTVar (const std::vector< T > &kappas, const T &val, const Int_t &bin, const Int_t whichDirection) const =0 |
|
virtual Bool_t | testRangeValidity (const T &val, const Int_t whichDirection) const =0 |
|
virtual void | cropValueForRange (T &val, const Int_t whichDirection) const =0 |
|
virtual T | interpolateFcn (Int_t code, const char *rangeName=0) const =0 |
|
virtual Double_t | evaluate () const =0 |
|
virtual void | getBArray (const std::vector< T > &kappas, const std::vector< T > &fcnList, std::vector< T > &BArray, BoundaryCondition const &bcBegin, BoundaryCondition const &bcEnd) const |
|
virtual void | getAArray (const std::vector< T > &kappas, std::vector< std::vector< T >> &AArray, BoundaryCondition const &bcBegin, BoundaryCondition const &bcEnd) const |
|
virtual std::vector< std::vector< T > > | getCoefficientsAlongDirection (const std::vector< T > &kappas, const TMatrix_t &Ainv, const std::vector< T > &fcnList, BoundaryCondition const &bcBegin, BoundaryCondition const &bcEnd, const Int_t pickBin) const |
|
virtual std::vector< T > | getCoefficients (const TVector_t &S, const std::vector< T > &kappas, const std::vector< T > &fcnList, const Int_t &bin) const |
|
virtual T | evalSplineSegment (const std::vector< T > &coefs, const T &kappa, const T &tup, const T &tdn, Bool_t doIntegrate=false) const |
|
Definition at line 16 of file MELANCSplineCore.h.
◆ TMatrix_t
◆ TVector_t
◆ BoundaryCondition
Enumerator |
---|
bcApproximatedSlope | |
bcClamped | |
bcApproximatedSecondDerivative | |
bcNaturalSpline | |
bcQuadratic | |
bcQuadraticWithNullSlope | |
NBoundaryConditions | |
Definition at line 28 of file MELANCSplineCore.h.
◆ VerbosityLevel
◆ MELANCSplineCore() [1/4]
MELANCSplineCore::MELANCSplineCore |
( |
| ) |
|
◆ MELANCSplineCore() [2/4]
MELANCSplineCore::MELANCSplineCore |
( |
const char * |
name, |
|
|
const char * |
title |
|
) |
| |
◆ MELANCSplineCore() [3/4]
MELANCSplineCore::MELANCSplineCore |
( |
const char * |
name, |
|
|
const char * |
title, |
|
|
RooAbsReal & |
inXVar, |
|
|
const std::vector< T > & |
inXList, |
|
|
Bool_t |
inUseFloor = true , |
|
|
T |
inFloorEval = 1e-15 , |
|
|
T |
inFloorInt = 1e-10 |
|
) |
| |
◆ MELANCSplineCore() [4/4]
MELANCSplineCore::MELANCSplineCore |
( |
const MELANCSplineCore & |
other, |
|
|
const char * |
name = 0 |
|
) |
| |
◆ ~MELANCSplineCore()
virtual MELANCSplineCore::~MELANCSplineCore |
( |
| ) |
|
|
inlinevirtual |
◆ addLeafDependents()
void MELANCSplineCore::addLeafDependents |
( |
RooArgSet & |
set | ) |
|
|
protected |
Definition at line 280 of file MELANCSplineCore.cc.
281 TIterator* iter =
set.createIterator();
283 while ((absarg = (RooAbsArg*)iter->Next())){
if (
dynamic_cast<RooRealVar*
>(absarg))
leafDepsList.add(*absarg); }
◆ analyticalIntegral()
virtual Double_t MELANCSplineCore::analyticalIntegral |
( |
Int_t |
code, |
|
|
const char * |
rangeName = 0 |
|
) |
| const |
|
pure virtual |
◆ clone()
virtual TObject* MELANCSplineCore::clone |
( |
const char * |
newname | ) |
const |
|
pure virtual |
◆ cropValueForRange()
virtual void MELANCSplineCore::cropValueForRange |
( |
T & |
val, |
|
|
const Int_t |
whichDirection |
|
) |
| const |
|
protectedpure virtual |
◆ doFloor()
void MELANCSplineCore::doFloor |
( |
Bool_t |
flag | ) |
|
◆ emptyFcnList()
virtual void MELANCSplineCore::emptyFcnList |
( |
| ) |
|
|
protectedpure virtual |
◆ evalSplineSegment()
MELANCSplineCore::T MELANCSplineCore::evalSplineSegment |
( |
const std::vector< T > & |
coefs, |
|
|
const T & |
kappa, |
|
|
const T & |
tup, |
|
|
const T & |
tdn, |
|
|
Bool_t |
doIntegrate = false |
|
) |
| const |
|
protectedvirtual |
Definition at line 265 of file MELANCSplineCore.cc.
267 for (
unsigned int ic=0; ic<coefs.size(); ic++){
268 if (doIntegrate) res += coefs.at(ic)*(pow(tup, (
int)(ic+1))-pow(tdn, (
int)(ic+1)))/((
MELANCSplineCore::T)(ic+1));
269 else res += coefs.at(ic)*pow(tup, (
int)ic);
271 if (doIntegrate) res /=
kappa;
◆ evaluate()
virtual Double_t MELANCSplineCore::evaluate |
( |
| ) |
const |
|
protectedpure virtual |
◆ getAArray()
void MELANCSplineCore::getAArray |
( |
const std::vector< T > & |
kappas, |
|
|
std::vector< std::vector< T >> & |
AArray, |
|
|
BoundaryCondition const & |
bcBegin, |
|
|
BoundaryCondition const & |
bcEnd |
|
) |
| const |
|
protectedvirtual |
Definition at line 159 of file MELANCSplineCore.cc.
161 Int_t npoints = kappas.size();
162 for (
int i=0; i<npoints; i++){
163 vector<MELANCSplineCore::T> Ai(npoints, 0.);
164 if (npoints==1) Ai[0]=1;
169 Ai[0]=2; Ai[1]=kappas.at(1)/kappas.at(0);
177 Ai[0]=1; Ai[1]=kappas.at(1)/kappas.at(0);
180 cerr <<
"MELANCSplineCore::getAArray: bcBegin " << bcBegin <<
" is not implemented!" << endl;
184 else if (i==npoints-1){
188 Ai[npoints-2]=1; Ai[npoints-1]=2.*kappas.at(npoints-1)/kappas.at(npoints-2);
196 Ai[npoints-2]=1; Ai[npoints-1]=kappas.at(npoints-1)/kappas.at(npoints-2);
199 cerr <<
"MELANCSplineCore::getAArray: bcEnd " << bcEnd <<
" is not implemented!" << endl;
211 Ai[i]=2.*kappa_j/kappa_jmo*(1.+kappa_j/kappa_jmo);
212 Ai[i+1]=kappa_j*kappa_jpo/pow(kappa_jmo, 2);
215 AArray.push_back(Ai);
◆ getAnalyticalIntegral()
virtual Int_t MELANCSplineCore::getAnalyticalIntegral |
( |
RooArgSet & |
allVars, |
|
|
RooArgSet & |
analVars, |
|
|
const char * |
rangeName = 0 |
|
) |
| const |
|
pure virtual |
◆ getBArray()
void MELANCSplineCore::getBArray |
( |
const std::vector< T > & |
kappas, |
|
|
const std::vector< T > & |
fcnList, |
|
|
std::vector< T > & |
BArray, |
|
|
BoundaryCondition const & |
bcBegin, |
|
|
BoundaryCondition const & |
bcEnd |
|
) |
| const |
|
protectedvirtual |
Definition at line 72 of file MELANCSplineCore.cc.
74 int npoints=kappas.size();
75 if (npoints!=(
int)fcnList.size()){
76 coutE(InputArguments) <<
"MELANCSplineCore::getBArray: Dim(kappas)=" << npoints <<
" != Dim(fcnList)=" << fcnList.size() << endl;
86 cerr <<
"Npoints " << npoints <<
" <= " << nthr << endl;
95 bcval=(npoints<3 ? 0 : ((fcnList.at(2)-fcnList.at(1))*kappas.at(1)-(fcnList.at(1)-fcnList.at(0))*kappas.at(0))/(0.5/kappas.at(0)+0.5/kappas.at(1)));
97 BArray.push_back(3.*(fcnList.at(1)-fcnList.at(0)) - bcval/2./pow(kappas.at(0), 2));
100 bcval=(fcnList.at(1)-fcnList.at(0))*kappas.at(0);
102 BArray.push_back(bcval/kappas.at(0));
105 bcval=2.*(fcnList.at(1)-fcnList.at(0));
106 BArray.push_back(bcval);
109 bcval=2.*(fcnList.at(1)-fcnList.at(0));
111 BArray.push_back(bcval*kappas.at(0)/kappas.at(1));
114 cerr <<
"MELANCSplineCore::getBArray: bcBegin " << bcBegin <<
" is not implemented!" << endl;
118 for (
int j=1; j<npoints-1; j++){
129 BArray.push_back(Bval);
134 ecval=(npoints<3 ? 0 : ((fcnList.at(npoints-1)-fcnList.at(npoints-2))*kappas.at(npoints-2)-(fcnList.at(npoints-2)-fcnList.at(npoints-3))*kappas.at(npoints-3))/(0.5/kappas.at(npoints-3)+0.5/kappas.at(npoints-2)));
136 BArray.push_back(3.*(fcnList.at(npoints-1)-fcnList.at(npoints-2)) + ecval/2./pow(kappas.at(npoints-1), 2));
139 ecval=(fcnList.at(npoints-1)-fcnList.at(npoints-2))*kappas.at(npoints-2);
141 BArray.push_back(ecval/kappas.at(npoints-1));
144 ecval=2.*(fcnList.at(npoints-1)-fcnList.at(npoints-2));
145 BArray.push_back(ecval);
148 ecval=2.*(fcnList.at(npoints-1)-fcnList.at(npoints-2));
149 BArray.push_back(ecval);
153 cerr <<
"MELANCSplineCore::getAArray: bcEnd " << bcEnd <<
" is not implemented!" << endl;
157 else if (npoints==1) BArray.push_back(0);
◆ getCoefficients()
vector< MELANCSplineCore::T > MELANCSplineCore::getCoefficients |
( |
const TVector_t & |
S, |
|
|
const std::vector< T > & |
kappas, |
|
|
const std::vector< T > & |
fcnList, |
|
|
const Int_t & |
bin |
|
) |
| const |
|
protectedvirtual |
Definition at line 219 of file MELANCSplineCore.cc.
221 vector<MELANCSplineCore::T> res;
223 const int fcnsize = fcnList.size();
227 if (fcnsize>(bin+1)){
233 C -= S[bin+1]*kappas.at(bin+1)/kappas.at(bin);
237 D += S[bin+1]*kappas.at(bin+1)/kappas.at(bin);
◆ getCoefficientsAlongDirection()
Definition at line 247 of file MELANCSplineCore.cc.
248 vector<MELANCSplineCore::T> BArray;
249 getBArray(kappas, fcnList, BArray, bcBegin, bcEnd);
251 Int_t npoints = BArray.size();
253 for (
int i=0; i<npoints; i++) Btrans[i]=BArray.at(i);
256 vector<vector<MELANCSplineCore::T>> coefs;
257 for (Int_t bin=0; bin<(npoints>1 ? npoints-1 : 1); bin++){
258 if (pickBin>=0 && bin!=pickBin)
continue;
259 vector<MELANCSplineCore::T> coef =
getCoefficients(Strans, kappas, fcnList, bin);
260 coefs.push_back(coef);
◆ getKappas()
virtual void MELANCSplineCore::getKappas |
( |
std::vector< T > & |
kappas, |
|
|
const Int_t |
whichDirection |
|
) |
| |
|
protectedpure virtual |
◆ getLeafDependents()
void MELANCSplineCore::getLeafDependents |
( |
RooRealProxy & |
proxy, |
|
|
RooArgSet & |
set |
|
) |
| |
|
protected |
◆ getTVar()
virtual T MELANCSplineCore::getTVar |
( |
const std::vector< T > & |
kappas, |
|
|
const T & |
val, |
|
|
const Int_t & |
bin, |
|
|
const Int_t |
whichDirection |
|
) |
| const |
|
protectedpure virtual |
◆ getWhichBin()
virtual Int_t MELANCSplineCore::getWhichBin |
( |
const T & |
val, |
|
|
const Int_t |
whichDirection |
|
) |
| const |
|
protectedpure virtual |
◆ interpolateFcn()
virtual T MELANCSplineCore::interpolateFcn |
( |
Int_t |
code, |
|
|
const char * |
rangeName = 0 |
|
) |
| const |
|
protectedpure virtual |
◆ npointsX()
unsigned int MELANCSplineCore::npointsX |
( |
| ) |
const |
|
inlineprotected |
◆ setEvalFloor()
◆ setIntFloor()
◆ setRangeValidity()
virtual void MELANCSplineCore::setRangeValidity |
( |
const T |
valmin, |
|
|
const T |
valmax, |
|
|
const Int_t |
whichDirection |
|
) |
| |
|
pure virtual |
◆ setVerbosity()
◆ testRangeValidity()
virtual Bool_t MELANCSplineCore::testRangeValidity |
( |
const T & |
val, |
|
|
const Int_t |
whichDirection |
|
) |
| const |
|
protectedpure virtual |
◆ floorEval
T MELANCSplineCore::floorEval |
|
protected |
◆ floorInt
T MELANCSplineCore::floorInt |
|
protected |
◆ leafDepsList
RooListProxy MELANCSplineCore::leafDepsList |
|
protected |
◆ rangeXmax
T MELANCSplineCore::rangeXmax |
|
protected |
◆ rangeXmin
T MELANCSplineCore::rangeXmin |
|
protected |
◆ theXVar
RooRealProxy MELANCSplineCore::theXVar |
|
protected |
◆ useFloor
Bool_t MELANCSplineCore::useFloor |
|
protected |
◆ verbosity
◆ XList
std::vector<T> MELANCSplineCore::XList |
|
protected |
The documentation for this class was generated from the following files: