#include "Riostream.h"
#include "RooParamHistFunc.h"
#include "RooAbsReal.h"
#include "RooAbsCategory.h"
#include "RooRealVar.h"
#include <math.h>
#include "TMath.h"
using namespace std ;
ClassImp(RooParamHistFunc)
;
RooParamHistFunc::RooParamHistFunc(const char *name, const char *title, RooDataHist& dh, Bool_t paramRelative) :
RooAbsReal(name,title),
_x("x","x",this),
_p("p","p",this),
_dh(dh),
_relParam(paramRelative)
{
_x.add(*_dh.get()) ;
RooArgSet allVars ;
for (Int_t i=0 ; i<_dh.numEntries() ; i++) {
_dh.get(i) ;
const char* vname = Form("%s_gamma_bin_%i",GetName(),i) ;
RooRealVar* var = new RooRealVar(vname,vname,0,1000) ;
var->setVal(_relParam ? 1 : _dh.weight()) ;
var->setConstant(kTRUE) ;
allVars.add(*var) ;
_p.add(*var) ;
}
addOwnedComponents(allVars) ;
}
RooParamHistFunc::RooParamHistFunc(const char *name, const char *title, const RooAbsArg& , RooDataHist& dh, Bool_t paramRelative) :
RooAbsReal(name,title),
_x("x","x",this),
_p("p","p",this),
_dh(dh),
_relParam(paramRelative)
{
_x.add(*_dh.get()) ;
RooArgSet allVars ;
for (Int_t i=0 ; i<_dh.numEntries() ; i++) {
_dh.get(i) ;
const char* vname = Form("%s_gamma_bin_%i",GetName(),i) ;
RooRealVar* var = new RooRealVar(vname,vname,0,1000) ;
var->setVal(_relParam ? 1 : _dh.weight()) ;
var->setConstant(kTRUE) ;
allVars.add(*var) ;
_p.add(*var) ;
}
addOwnedComponents(allVars) ;
}
RooParamHistFunc::RooParamHistFunc(const char *name, const char *title, RooDataHist& dh, const RooParamHistFunc& paramSource, Bool_t paramRelative) :
RooAbsReal(name,title),
_x("x","x",this),
_p("p","p",this),
_dh(dh),
_relParam(paramRelative)
{
_x.add(*_dh.get()) ;
_p.add(paramSource._p) ;
}
RooParamHistFunc::RooParamHistFunc(const RooParamHistFunc& other, const char* name) :
RooAbsReal(other,name),
_x("x",this,other._x),
_p("p",this,other._p),
_dh(other._dh),
_relParam(other._relParam)
{
}
Double_t RooParamHistFunc::evaluate() const
{
Int_t idx = ((RooDataHist&)_dh).getIndex(_x,kTRUE) ;
Double_t ret = ((RooAbsReal*)_p.at(idx))->getVal() ;
if (_relParam) {
Double_t nom = getNominal(idx) ;
ret *= nom ;
}
return ret ;
}
Double_t RooParamHistFunc::getActual(Int_t ibin)
{
return ((RooAbsReal&)_p[ibin]).getVal() ;
}
void RooParamHistFunc::setActual(Int_t ibin, Double_t newVal)
{
((RooRealVar&)_p[ibin]).setVal(newVal) ;
}
Double_t RooParamHistFunc::getNominal(Int_t ibin) const
{
_dh.get(ibin) ;
return _dh.weight() ;
}
Double_t RooParamHistFunc::getNominalError(Int_t ibin) const
{
_dh.get(ibin) ;
return _dh.weightError() ;
}
list<Double_t>* RooParamHistFunc::plotSamplingHint(RooAbsRealLValue& obs, Double_t xlo, Double_t xhi) const
{
RooAbsLValue* lvarg = dynamic_cast<RooAbsLValue*>(_dh.get()->find(obs.GetName())) ;
if (!lvarg) {
return 0 ;
}
const RooAbsBinning* binning = lvarg->getBinningPtr(0) ;
Double_t* boundaries = binning->array() ;
list<Double_t>* hint = new list<Double_t> ;
xlo = xlo - 0.01*(xhi-xlo) ;
xhi = xhi + 0.01*(xhi-xlo) ;
Double_t delta = (xhi-xlo)*1e-8 ;
for (Int_t i=0 ; i<binning->numBoundaries() ; i++) {
if (boundaries[i]>=xlo && boundaries[i]<=xhi) {
hint->push_back(boundaries[i]-delta) ;
hint->push_back(boundaries[i]+delta) ;
}
}
return hint ;
}
std::list<Double_t>* RooParamHistFunc::binBoundaries(RooAbsRealLValue& obs, Double_t xlo, Double_t xhi) const
{
RooAbsLValue* lvarg = dynamic_cast<RooAbsLValue*>(_dh.get()->find(obs.GetName())) ;
if (!lvarg) {
return 0 ;
}
const RooAbsBinning* binning = lvarg->getBinningPtr(0) ;
Double_t* boundaries = binning->array() ;
list<Double_t>* hint = new list<Double_t> ;
for (Int_t i=0 ; i<binning->numBoundaries() ; i++) {
if (boundaries[i]>=xlo && boundaries[i]<=xhi) {
hint->push_back(boundaries[i]) ;
}
}
return hint ;
}
Int_t RooParamHistFunc::getAnalyticalIntegralWN(RooArgSet& allVars, RooArgSet& analVars,
const RooArgSet* , const char* ) const
{
RooAbsCollection *allVarsCommon = allVars.selectCommon(_x) ;
Bool_t intAllObs = (allVarsCommon->getSize()==_x.getSize()) ;
delete allVarsCommon ;
if (intAllObs && matchArgs(allVars,analVars,_x)) {
return 1 ;
}
return 0 ;
}
Double_t RooParamHistFunc::analyticalIntegralWN(Int_t code, const RooArgSet* ,const char* ) const
{
R__ASSERT(code==1) ;
RooFIter iter = _p.fwdIterator() ;
RooAbsReal* p ;
Double_t ret(0) ;
Int_t i(0) ;
while((p=(RooAbsReal*)iter.next())) {
Double_t bin = p->getVal() ;
if (_relParam) bin *= getNominal(i++) ;
ret += bin ;
}
RooFIter xiter = _x.fwdIterator() ;
RooAbsArg* obs ;
Double_t binV(1) ;
while((obs=xiter.next())) {
RooRealVar* xx = (RooRealVar*) obs ;
binV *= (xx->getMax()-xx->getMin())/xx->numBins() ;
}
return ret*binV ;
}