#include "RooFit.h"
#include "Riostream.h"
#include "Riostream.h"
#include <math.h>
#include "RooParametricStepFunction.h"
#include "RooAbsReal.h"
#include "RooRealVar.h"
#include "RooArgList.h"
ClassImp(RooParametricStepFunction)
;
RooParametricStepFunction::RooParametricStepFunction(const char* name, const char* title,
RooAbsReal& x, const RooArgList& coefList, TArrayD& limits, Int_t nBins) :
RooAbsPdf(name, title),
_x("x", "Dependent", this, x),
_coefList("coefList","List of coefficients",this),
_nBins(nBins)
{
_coefIter = _coefList.createIterator() ;
if (_nBins<0) {
cout << "RooParametricStepFunction::ctor(" << GetName()
<< ") WARNING: nBins must be >=0, setting value to 0" << endl ;
_nBins=0 ;
}
TIterator* coefIter = coefList.createIterator() ;
RooAbsArg* coef ;
while((coef = (RooAbsArg*)coefIter->Next())) {
if (!dynamic_cast<RooAbsReal*>(coef)) {
cout << "RooParametricStepFunction::ctor(" << GetName() << ") ERROR: coefficient " << coef->GetName()
<< " is not of type RooAbsReal" << endl ;
assert(0) ;
}
_coefList.add(*coef) ;
}
delete coefIter ;
limits.Copy(_limits);
}
RooParametricStepFunction::RooParametricStepFunction(const RooParametricStepFunction& other, const char* name) :
RooAbsPdf(other, name),
_x("x", this, other._x),
_coefList("coefList",this,other._coefList),
_nBins(other._nBins)
{
_coefIter = _coefList.createIterator();
(other._limits).Copy(_limits);
}
RooParametricStepFunction::~RooParametricStepFunction()
{
delete _coefIter ;
}
Int_t RooParametricStepFunction::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* ) const
{
if (matchArgs(allVars, analVars, _x)) return 1;
return 0;
}
Double_t RooParametricStepFunction::analyticalIntegral(Int_t code, const char* rangeName) const
{
assert(code==1) ;
if (!rangeName) {
return 1.0 ;
}
Double_t xmin = _x.min(rangeName) ;
Double_t xmax = _x.max(rangeName) ;
Double_t sum=0 ;
Int_t i ;
for (i=1 ; i<=_nBins ; i++) {
Double_t binVal = (i<_nBins) ? (static_cast<RooRealVar*>(_coefList.at(i-1))->getVal()) : lastBinValue() ;
if (_limits[i-1]>=xmin && _limits[i]<=xmax) {
sum += (_limits[i]-_limits[i-1])*binVal ;
} else if (_limits[i-1]<xmin && _limits[i]>xmax) {
sum += (xmax-xmin)*binVal ;
return sum ;
} else if (_limits[i-1]<xmin && _limits[i]<=xmax && _limits[i]>xmin) {
sum += (_limits[i]-xmin)*binVal ;
} else if (_limits[i-1]>=xmin && _limits[i]>xmax && _limits[i-1]<xmax) {
sum += (xmax-_limits[i-1])*binVal ;
return sum ;
}
}
return sum;
}
Double_t RooParametricStepFunction::lastBinValue() const
{
Double_t sum(0.);
Double_t binSize(0.);
for (Int_t j=1;j<_nBins;j++){
RooRealVar* tmp = (RooRealVar*) _coefList.at(j-1);
binSize = _limits[j] - _limits[j-1];
sum = sum + tmp->getVal()*binSize;
}
binSize = _limits[_nBins] - _limits[_nBins-1];
return (1.0 - sum)/binSize;
}
Double_t RooParametricStepFunction::evaluate() const
{
Double_t xval(0.);
xval = _x;
Double_t value(0.);
if (_x >= _limits[0] && _x < _limits[_nBins]){
for (Int_t i=1;i<=_nBins;i++){
if (_x < _limits[i]){
if (i<_nBins) {
RooRealVar* tmp = (RooRealVar*) _coefList.at(i-1);
value = tmp->getVal();
break;
} else {
Double_t sum(0.);
Double_t binSize(0.);
for (Int_t j=1;j<_nBins;j++){
RooRealVar* tmp = (RooRealVar*) _coefList.at(j-1);
binSize = _limits[j] - _limits[j-1];
sum = sum + tmp->getVal()*binSize;
}
binSize = _limits[_nBins] - _limits[_nBins-1];
value = (1.0 - sum)/binSize;
if (value<=0.0){
value = 0.000001;
}
break;
}
}
}
}
return value;
}
Int_t RooParametricStepFunction::getnBins(){
return _nBins;
}
Double_t* RooParametricStepFunction::getLimits(){
Double_t* limoutput = _limits.GetArray();
return limoutput;
}
ROOT page - Class index - Class Hierarchy - Top of the page
This page has been automatically generated. If you have any comments or suggestions about the page layout send a mail to ROOT support, or contact the developers with any questions or problems regarding ROOT.