#include "RooFit.h"
#include "TClass.h"
#include "TMath.h"
#include "Riostream.h"
#include "RooResolutionModel.h"
ClassImp(RooResolutionModel)
;
RooFormulaVar* RooResolutionModel::_identity = 0;
RooResolutionModel::RooResolutionModel(const char *name, const char *title, RooRealVar& _x) :
RooAbsPdf(name,title),
x("x","Dependent or convolution variable",this,_x),
_basisCode(0), _basis(0),
_ownBasis(kFALSE)
{
if (!_identity) _identity = new RooFormulaVar("identity","1",RooArgSet("")) ;
}
RooResolutionModel::RooResolutionModel(const RooResolutionModel& other, const char* name) :
RooAbsPdf(other,name),
x("x",this,other.x),
_basisCode(other._basisCode), _basis(0),
_ownBasis(kFALSE)
{
if (other._basis) {
_basis = (RooFormulaVar*) other._basis->Clone() ;
_ownBasis = kTRUE ;
}
if (_basis) {
TIterator* bsIter = _basis->serverIterator() ;
RooAbsArg* basisServer ;
while((basisServer = (RooAbsArg*)bsIter->Next())) {
addServer(*basisServer,kTRUE,kFALSE) ;
}
delete bsIter ;
}
}
RooResolutionModel::~RooResolutionModel()
{
if (_ownBasis && _basis) {
delete _basis ;
}
}
RooFormulaVar* RooResolutionModel::identity()
{
return _identity ;
}
RooResolutionModel* RooResolutionModel::convolution(RooFormulaVar* basis, RooAbsArg* owner) const
{
if (basis->getParameter(0) != x.absArg()) {
cout << "RooResolutionModel::convolution(" << GetName() << "," << this
<< ") convolution parameter of basis function and PDF don't match" << endl ;
cout << "basis->findServer(0) = " << basis->findServer(0) << endl ;
cout << "x.absArg() = " << x.absArg() << endl ;
return 0 ;
}
TString newName(GetName()) ;
newName.Append("_conv_") ;
newName.Append(basis->GetName()) ;
newName.Append("_[") ;
newName.Append(owner->GetName()) ;
newName.Append("]") ;
RooResolutionModel* conv = (RooResolutionModel*) clone(newName) ;
TString newTitle(conv->GetTitle()) ;
newTitle.Append(" convoluted with basis function ") ;
newTitle.Append(basis->GetName()) ;
conv->SetTitle(newTitle.Data()) ;
conv->changeBasis(basis) ;
return conv ;
}
void RooResolutionModel::changeBasis(RooFormulaVar* basis)
{
if (_basis) {
TIterator* bsIter = _basis->serverIterator() ;
RooAbsArg* basisServer ;
while((basisServer = (RooAbsArg*)bsIter->Next())) {
removeServer(*basisServer) ;
}
delete bsIter ;
if (_ownBasis) {
delete _basis ;
}
}
_ownBasis = kFALSE ;
_basis = basis ;
if (_basis) {
TIterator* bsIter = _basis->serverIterator() ;
RooAbsArg* basisServer ;
while((basisServer = (RooAbsArg*)bsIter->Next())) {
addServer(*basisServer,kTRUE,kFALSE) ;
}
delete bsIter ;
}
_basisCode = basis?basisCode(basis->GetTitle()):0 ;
}
const RooRealVar& RooResolutionModel::basisConvVar() const
{
TIterator* sIter = basis().serverIterator() ;
RooRealVar* var = (RooRealVar*) sIter->Next() ;
delete sIter ;
return *var ;
}
RooRealVar& RooResolutionModel::convVar() const
{
return (RooRealVar&) x.arg() ;
}
Double_t RooResolutionModel::getVal(const RooArgSet* nset) const
{
if (!_basis) return RooAbsPdf::getVal(nset) ;
if (isValueDirty()) {
_value = evaluate() ;
if (_verboseDirty) cout << "RooResolutionModel(" << GetName() << ") value = " << _value << endl ;
clearValueDirty() ;
clearShapeDirty() ;
}
return _value ;
}
Bool_t RooResolutionModel::redirectServersHook(const RooAbsCollection& newServerList, Bool_t mustReplaceAll, Bool_t nameChange, Bool_t )
{
if (!_basis) return kFALSE ;
RooFormulaVar* newBasis = (RooFormulaVar*) newServerList.find(_basis->GetName()) ;
if (newBasis) {
if (_ownBasis) {
delete _basis ;
}
_basis = newBasis ;
_ownBasis = kFALSE ;
}
_basis->redirectServers(newServerList,mustReplaceAll,nameChange) ;
return (mustReplaceAll && !newBasis) ;
}
Bool_t RooResolutionModel::traceEvalHook(Double_t value) const
{
return isnan(value) ;
}
void RooResolutionModel::normLeafServerList(RooArgSet& list) const
{
_norm->leafNodeServerList(&list) ;
}
Double_t RooResolutionModel::getNorm(const RooArgSet* nset) const
{
if (!nset) {
return getVal() ;
}
syncNormalization(nset,kFALSE) ;
if (_verboseEval>1) cout << IsA()->GetName() << "::getNorm(" << GetName()
<< "): norm(" << _norm << ") = " << _norm->getVal() << endl ;
Double_t ret = _norm->getVal() ;
return ret ;
}
void RooResolutionModel::printToStream(ostream& os, PrintOption opt, TString indent) const
{
RooAbsPdf::printToStream(os,opt,indent) ;
if(opt >= Verbose) {
os << indent << "--- RooResolutionModel ---" << endl;
os << indent << "basis function = " ;
if (_basis) {
_basis->printToStream(os,opt,indent) ;
} else {
os << "<none>" << endl ;
}
}
}
Bool_t RooResolutionModel::syncNormalizationPreHook(RooAbsReal*,const RooArgSet*) const
{
return (_basisCode!=0) ;
}
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.