#include "RooFit.h"
#include "RooMsgService.h"
#include "TClass.h"
#include "Riostream.h"
#include "TMath.h"
#include "TObjString.h"
#include "TPaveText.h"
#include "TList.h"
#include "TH1.h"
#include "TH2.h"
#include "TMatrixD.h"
#include "TMatrixDSym.h"
#include "RooAbsPdf.h"
#include "RooDataSet.h"
#include "RooArgSet.h"
#include "RooArgProxy.h"
#include "RooRealProxy.h"
#include "RooRealVar.h"
#include "RooGenContext.h"
#include "RooBinnedGenContext.h"
#include "RooPlot.h"
#include "RooCurve.h"
#include "RooNLLVar.h"
#include "RooMinuit.h"
#include "RooCategory.h"
#include "RooNameReg.h"
#include "RooCmdConfig.h"
#include "RooGlobalFunc.h"
#include "RooAddition.h"
#include "RooRandom.h"
#include "RooNumIntConfig.h"
#include "RooProjectedPdf.h"
#include "RooInt.h"
#include "RooCustomizer.h"
#include "RooConstraintSum.h"
#include "RooParamBinning.h"
#include "RooNumCdf.h"
#include "RooFitResult.h"
#include "RooNumGenConfig.h"
#include "RooCachedReal.h"
#include "RooXYChi2Var.h"
#include "RooChi2Var.h"
#include "RooMinimizer.h"
#include "RooRealIntegral.h"
#include "Math/CholeskyDecomp.h"
#include <string>
using namespace std;
ClassImp(RooAbsPdf)
;
ClassImp(RooAbsPdf::GenSpec)
;
Int_t RooAbsPdf::_verboseEval = 0;
Bool_t RooAbsPdf::_evalError = kFALSE ;
TString RooAbsPdf::_normRangeOverride ;
RooAbsPdf::RooAbsPdf() : _norm(0), _normSet(0), _specGeneratorConfig(0)
{
_errorCount = 0 ;
_negCount = 0 ;
_rawValue = 0 ;
_selectComp = kFALSE ;
_traceCount = 0 ;
}
RooAbsPdf::RooAbsPdf(const char *name, const char *title) :
RooAbsReal(name,title), _norm(0), _normSet(0), _normMgr(this,10), _selectComp(kTRUE), _specGeneratorConfig(0)
{
resetErrorCounters() ;
setTraceCounter(0) ;
}
RooAbsPdf::RooAbsPdf(const char *name, const char *title,
Double_t plotMin, Double_t plotMax) :
RooAbsReal(name,title,plotMin,plotMax), _norm(0), _normSet(0), _normMgr(this,10), _selectComp(kTRUE), _specGeneratorConfig(0)
{
resetErrorCounters() ;
setTraceCounter(0) ;
}
RooAbsPdf::RooAbsPdf(const RooAbsPdf& other, const char* name) :
RooAbsReal(other,name), _norm(0), _normSet(0),
_normMgr(other._normMgr,this), _selectComp(other._selectComp), _normRange(other._normRange)
{
resetErrorCounters() ;
setTraceCounter(other._traceCount) ;
if (other._specGeneratorConfig) {
_specGeneratorConfig = new RooNumGenConfig(*other._specGeneratorConfig) ;
} else {
_specGeneratorConfig = 0 ;
}
}
RooAbsPdf::~RooAbsPdf()
{
if (_specGeneratorConfig) delete _specGeneratorConfig ;
}
Double_t RooAbsPdf::getValV(const RooArgSet* nset) const
{
if (!nset) {
RooArgSet* tmp = _normSet ;
_normSet = 0 ;
Double_t val = evaluate() ;
_normSet = tmp ;
Bool_t error = traceEvalPdf(val) ;
if (error) {
return 0 ;
}
return val ;
}
Bool_t nsetChanged(kFALSE) ;
if (nset!=_normSet || _norm==0) {
nsetChanged = syncNormalization(nset) ;
}
if (isValueDirty() || nsetChanged || _norm->isValueDirty()) {
Double_t rawVal = evaluate() ;
Bool_t error = traceEvalPdf(rawVal) ;
Double_t normVal(_norm->getVal()) ;
if (normVal<=0.) {
error=kTRUE ;
logEvalError("p.d.f normalization integral is zero or negative") ;
}
if (error) {
_value = 0 ;
} else {
_value = rawVal / normVal ;
}
clearValueAndShapeDirty() ;
}
return _value ;
}
Double_t RooAbsPdf::analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* rangeName) const
{
cxcoutD(Eval) << "RooAbsPdf::analyticalIntegralWN(" << GetName() << ") code = " << code << " normset = " << (normSet?*normSet:RooArgSet()) << endl ;
if (code==0) return getVal(normSet) ;
if (normSet) {
return analyticalIntegral(code,rangeName) / getNorm(normSet) ;
} else {
return analyticalIntegral(code,rangeName) ;
}
}
Bool_t RooAbsPdf::traceEvalPdf(Double_t value) const
{
Bool_t error(kFALSE) ;
if (TMath::IsNaN(value)) {
logEvalError(Form("p.d.f value is Not-a-Number (%f), forcing value to zero",value)) ;
error=kTRUE ;
}
if (value<0) {
logEvalError(Form("p.d.f value is less than zero (%f), forcing value to zero",value)) ;
error=kTRUE ;
}
if(!error) return error ;
if(++_errorCount <= 10) {
cxcoutD(Tracing) << "*** Evaluation Error " << _errorCount << " ";
if(_errorCount == 10) cxcoutD(Tracing) << "(no more will be printed) ";
}
else {
return error ;
}
Print() ;
return error ;
}
Double_t RooAbsPdf::getNorm(const RooArgSet* nset) const
{
if (!nset) return 1 ;
syncNormalization(nset,kTRUE) ;
if (_verboseEval>1) cxcoutD(Tracing) << IsA()->GetName() << "::getNorm(" << GetName() << "): norm(" << _norm << ") = " << _norm->getVal() << endl ;
Double_t ret = _norm->getVal() ;
if (ret==0.) {
if(++_errorCount <= 10) {
coutW(Eval) << "RooAbsPdf::getNorm(" << GetName() << ":: WARNING normalization is zero, nset = " ; nset->Print("1") ;
if(_errorCount == 10) coutW(Eval) << "RooAbsPdf::getNorm(" << GetName() << ") INFO: no more messages will be printed " << endl ;
}
}
return ret ;
}
const RooAbsReal* RooAbsPdf::getNormObj(const RooArgSet* nset, const RooArgSet* iset, const TNamed* rangeName) const
{
CacheElem* cache = (CacheElem*) _normMgr.getObj(nset,iset,0,rangeName) ;
if (cache) {
return cache->_norm ;
}
RooArgSet* depList = getObservables(iset) ;
RooAbsReal* norm = createIntegral(*depList,*nset, *getIntegratorConfig(), RooNameReg::str(rangeName)) ;
delete depList ;
cache = new CacheElem(*norm) ;
_normMgr.setObj(nset,iset,cache,rangeName) ;
return norm ;
}
Bool_t RooAbsPdf::syncNormalization(const RooArgSet* nset, Bool_t adjustProxies) const
{
_normSet = (RooArgSet*) nset ;
CacheElem* cache = (CacheElem*) _normMgr.getObj(nset) ;
if (cache) {
Bool_t nsetChanged = (_norm!=cache->_norm) ;
_norm = cache->_norm ;
if (nsetChanged && adjustProxies) {
((RooAbsPdf*) this)->setProxyNormSet(nset) ;
}
return nsetChanged ;
}
if (adjustProxies) {
((RooAbsPdf*) this)->setProxyNormSet(nset) ;
}
RooArgSet* depList = getObservables(nset) ;
if (_verboseEval>0) {
if (!selfNormalized()) {
cxcoutD(Tracing) << IsA()->GetName() << "::syncNormalization(" << GetName()
<< ") recreating normalization integral " << endl ;
if (depList) depList->printStream(ccoutD(Tracing),kName|kValue|kArgs,kSingleLine) ; else ccoutD(Tracing) << "<none>" << endl ;
} else {
cxcoutD(Tracing) << IsA()->GetName() << "::syncNormalization(" << GetName() << ") selfNormalized, creating unit norm" << endl;
}
}
if (selfNormalized() || !dependsOn(*depList)) {
TString ntitle(GetTitle()) ; ntitle.Append(" Unit Normalization") ;
TString nname(GetName()) ; nname.Append("_UnitNorm") ;
_norm = new RooRealVar(nname.Data(),ntitle.Data(),1) ;
} else {
const char* nr = (_normRangeOverride.Length()>0 ? _normRangeOverride.Data() : (_normRange.Length()>0 ? _normRange.Data() : 0)) ;
RooAbsReal* normInt = createIntegral(*depList,*getIntegratorConfig(),nr) ;
normInt->getVal() ;
const char* cacheParamsStr = getStringAttribute("CACHEPARAMINT") ;
if (cacheParamsStr && strlen(cacheParamsStr)) {
RooArgSet* intParams = normInt->getVariables() ;
RooNameSet cacheParamNames ;
cacheParamNames.setNameList(cacheParamsStr) ;
RooArgSet* cacheParams = cacheParamNames.select(*intParams) ;
if (cacheParams->getSize()>0) {
cxcoutD(Caching) << "RooAbsReal::createIntObj(" << GetName() << ") INFO: constructing " << cacheParams->getSize()
<< "-dim value cache for integral over " << *depList << " as a function of " << *cacheParams << " in range " << (nr?nr:"<default>") << endl ;
string name = Form("%s_CACHE_[%s]",normInt->GetName(),cacheParams->contentsString().c_str()) ;
RooCachedReal* cachedIntegral = new RooCachedReal(name.c_str(),name.c_str(),*normInt,*cacheParams) ;
cachedIntegral->setInterpolationOrder(2) ;
cachedIntegral->addOwnedComponents(*normInt) ;
cachedIntegral->setCacheSource(kTRUE) ;
if (normInt->operMode()==ADirty) {
cachedIntegral->setOperMode(ADirty) ;
}
normInt= cachedIntegral ;
}
delete cacheParams ;
delete intParams ;
}
_norm = normInt ;
}
cache = new CacheElem(*_norm) ;
_normMgr.setObj(nset,cache) ;
delete depList ;
return kTRUE ;
}
Bool_t RooAbsPdf::traceEvalHook(Double_t value) const
{
Bool_t error= TMath::IsNaN(value) || (value < 0);
if(!error && _traceCount <= 0) return error ;
if(error && ++_errorCount <= 10) {
cxcoutD(Tracing) << "*** Evaluation Error " << _errorCount << " ";
if(_errorCount == 10) ccoutD(Tracing) << "(no more will be printed) ";
}
else if(_traceCount > 0) {
ccoutD(Tracing) << '[' << _traceCount-- << "] ";
}
else {
return error ;
}
Print() ;
return error ;
}
void RooAbsPdf::resetErrorCounters(Int_t resetValue)
{
_errorCount = resetValue ;
_negCount = resetValue ;
}
void RooAbsPdf::setTraceCounter(Int_t value, Bool_t allNodes)
{
if (!allNodes) {
_traceCount = value ;
return ;
} else {
RooArgList branchList ;
branchNodeServerList(&branchList) ;
TIterator* iter = branchList.createIterator() ;
RooAbsArg* arg ;
while((arg=(RooAbsArg*)iter->Next())) {
RooAbsPdf* pdf = dynamic_cast<RooAbsPdf*>(arg) ;
if (pdf) pdf->setTraceCounter(value,kFALSE) ;
}
delete iter ;
}
}
Double_t RooAbsPdf::getLogVal(const RooArgSet* nset) const
{
Double_t prob = getVal(nset) ;
if (fabs(prob)>1e6) {
coutW(Eval) << "RooAbsPdf::getLogVal(" << GetName() << ") WARNING: large likelihood value: " << prob << endl ;
}
if(prob < 0) {
logEvalError("getLogVal() top-level p.d.f evaluates to a negative number") ;
return 0;
}
if(prob == 0) {
logEvalError("getLogVal() top-level p.d.f evaluates to zero") ;
return log((double)0);
}
if (TMath::IsNaN(prob)) {
logEvalError("getLogVal() top-level p.d.f evaluates to NaN") ;
return log((double)0);
}
return log(prob);
}
Double_t RooAbsPdf::extendedTerm(Double_t observed, const RooArgSet* nset) const
{
if(!canBeExtended()) {
coutE(InputArguments) << fName << ": this PDF does not support extended maximum likelihood"
<< endl;
return 0;
}
Double_t expected= expectedEvents(nset);
if(expected < 0) {
coutE(InputArguments) << fName << ": calculated negative expected events: " << expected
<< endl;
return 0;
}
if (fabs(expected)<1e-10 && fabs(observed)<1e-10) {
return 0 ;
}
if (expected<0 || TMath::IsNaN(expected)) {
logEvalError("extendedTerm #expected events is <0 or NaN") ;
return 0 ;
}
Double_t extra= expected - observed*log(expected);
Bool_t trace(kFALSE) ;
if(trace) {
cxcoutD(Tracing) << fName << "::extendedTerm: expected " << expected << " events, got "
<< observed << " events. extendedTerm = " << extra << endl;
}
return extra;
}
RooAbsReal* RooAbsPdf::createNLL(RooAbsData& data, const RooCmdArg& arg1, const RooCmdArg& arg2, const RooCmdArg& arg3, const RooCmdArg& arg4,
const RooCmdArg& arg5, const RooCmdArg& arg6, const RooCmdArg& arg7, const RooCmdArg& arg8)
{
RooLinkedList l ;
l.Add((TObject*)&arg1) ; l.Add((TObject*)&arg2) ;
l.Add((TObject*)&arg3) ; l.Add((TObject*)&arg4) ;
l.Add((TObject*)&arg5) ; l.Add((TObject*)&arg6) ;
l.Add((TObject*)&arg7) ; l.Add((TObject*)&arg8) ;
return createNLL(data,l) ;
}
RooAbsReal* RooAbsPdf::createNLL(RooAbsData& data, const RooLinkedList& cmdList)
{
RooCmdConfig pc(Form("RooAbsPdf::createNLL(%s)",GetName())) ;
pc.defineString("rangeName","RangeWithName",0,"",kTRUE) ;
pc.defineString("addCoefRange","SumCoefRange",0,"") ;
pc.defineString("globstag","GlobalObservablesTag",0,"") ;
pc.defineDouble("rangeLo","Range",0,-999.) ;
pc.defineDouble("rangeHi","Range",1,-999.) ;
pc.defineInt("splitRange","SplitRange",0,0) ;
pc.defineInt("ext","Extended",0,2) ;
pc.defineInt("numcpu","NumCPU",0,1) ;
pc.defineInt("interleave","NumCPU",1,0) ;
pc.defineInt("verbose","Verbose",0,0) ;
pc.defineInt("optConst","Optimize",0,0) ;
pc.defineInt("cloneData","CloneData",2,0) ;
pc.defineSet("projDepSet","ProjectedObservables",0,0) ;
pc.defineSet("cPars","Constrain",0,0) ;
pc.defineSet("glObs","GlobalObservables",0,0) ;
pc.defineInt("constrAll","Constrained",0,0) ;
pc.defineInt("doOffset","OffsetLikelihood",0,0) ;
pc.defineSet("extCons","ExternalConstraints",0,0) ;
pc.defineMutex("Range","RangeWithName") ;
pc.defineMutex("Constrain","Constrained") ;
pc.defineMutex("GlobalObservables","GlobalObservablesTag") ;
pc.process(cmdList) ;
if (!pc.ok(kTRUE)) {
return 0 ;
}
const char* rangeName = pc.getString("rangeName",0,kTRUE) ;
const char* addCoefRangeName = pc.getString("addCoefRange",0,kTRUE) ;
const char* globsTag = pc.getString("globstag",0,kTRUE) ;
Int_t ext = pc.getInt("ext") ;
Int_t numcpu = pc.getInt("numcpu") ;
RooFit::MPSplit interl = (RooFit::MPSplit) pc.getInt("interleave") ;
Int_t splitr = pc.getInt("splitRange") ;
Bool_t verbose = pc.getInt("verbose") ;
Int_t optConst = pc.getInt("optConst") ;
Int_t cloneData = pc.getInt("cloneData") ;
Int_t doOffset = pc.getInt("doOffset") ;
if (cloneData==2) {
cloneData = optConst ;
}
RooArgSet* cPars = pc.getSet("cPars") ;
RooArgSet* glObs = pc.getSet("glObs") ;
if (pc.hasProcessed("GlobalObservablesTag")) {
if (glObs) delete glObs ;
RooArgSet* allVars = getVariables() ;
glObs = (RooArgSet*) allVars->selectByAttrib(globsTag,kTRUE) ;
coutI(Minimization) << "User-defined specification of global observables definition with tag named '" << globsTag << "'" << endl ;
delete allVars ;
} else if (!pc.hasProcessed("GlobalObservables")) {
const char* defGlobObsTag = getStringAttribute("DefaultGlobalObservablesTag") ;
if (defGlobObsTag) {
coutI(Minimization) << "p.d.f. provides built-in specification of global observables definition with tag named '" << defGlobObsTag << "'" << endl ;
if (glObs) delete glObs ;
RooArgSet* allVars = getVariables() ;
glObs = (RooArgSet*) allVars->selectByAttrib(defGlobObsTag,kTRUE) ;
}
}
Bool_t doStripDisconnected=kFALSE ;
if (!cPars) {
cPars = getParameters(data,kFALSE) ;
doStripDisconnected=kTRUE ;
}
const RooArgSet* extCons = pc.getSet("extCons") ;
if (ext==2) {
ext = ((extendMode()==CanBeExtended || extendMode()==MustBeExtended)) ? 1 : 0 ;
if (ext) {
coutI(Minimization) << "p.d.f. provides expected number of events, including extended term in likelihood." << endl ;
}
}
if (pc.hasProcessed("Range")) {
Double_t rangeLo = pc.getDouble("rangeLo") ;
Double_t rangeHi = pc.getDouble("rangeHi") ;
RooArgSet* obs = getObservables(&data) ;
TIterator* iter = obs->createIterator() ;
RooAbsArg* arg ;
while((arg=(RooAbsArg*)iter->Next())) {
RooRealVar* rrv = dynamic_cast<RooRealVar*>(arg) ;
if (rrv) rrv->setRange("fit",rangeLo,rangeHi) ;
}
rangeName = "fit" ;
}
RooArgSet projDeps ;
RooArgSet* tmp = pc.getSet("projDepSet") ;
if (tmp) {
projDeps.add(*tmp) ;
}
RooAbsReal::setEvalErrorLoggingMode(RooAbsReal::CollectErrors) ;
RooAbsReal* nll ;
string baseName = Form("nll_%s_%s",GetName(),data.GetName()) ;
if (!rangeName || strchr(rangeName,',')==0) {
nll = new RooNLLVar(baseName.c_str(),"-log(likelihood)",*this,data,projDeps,ext,rangeName,addCoefRangeName,numcpu,interl,verbose,splitr,cloneData) ;
} else {
RooArgList nllList ;
const size_t bufSize = strlen(rangeName)+1;
char* buf = new char[bufSize] ;
strlcpy(buf,rangeName,bufSize) ;
char* token = strtok(buf,",") ;
while(token) {
RooAbsReal* nllComp = new RooNLLVar(Form("%s_%s",baseName.c_str(),token),"-log(likelihood)",*this,data,projDeps,ext,token,addCoefRangeName,numcpu,interl,verbose,splitr,cloneData) ;
nllList.add(*nllComp) ;
token = strtok(0,",") ;
}
delete[] buf ;
nll = new RooAddition(baseName.c_str(),"-log(likelihood)",nllList,kTRUE) ;
}
RooAbsReal::setEvalErrorLoggingMode(RooAbsReal::PrintErrors) ;
RooArgSet allConstraints ;
if (cPars && cPars->getSize()>0) {
RooArgSet* constraints = getAllConstraints(*data.get(),*cPars,doStripDisconnected) ;
allConstraints.add(*constraints) ;
delete constraints ;
}
if (extCons) {
allConstraints.add(*extCons) ;
}
RooAbsReal* nllCons(0) ;
if (allConstraints.getSize()>0 && cPars) {
coutI(Minimization) << " Including the following contraint terms in minimization: " << allConstraints << endl ;
if (glObs) {
coutI(Minimization) << "The following global observables have been defined: " << *glObs << endl ;
}
nllCons = new RooConstraintSum(Form("%s_constr",baseName.c_str()),"nllCons",allConstraints,glObs ? *glObs : *cPars) ;
nllCons->setOperMode(ADirty) ;
RooAbsReal* orignll = nll ;
nll = new RooAddition(Form("%s_with_constr",baseName.c_str()),"nllWithCons",RooArgSet(*nll,*nllCons)) ;
nll->addOwnedComponents(RooArgSet(*orignll,*nllCons)) ;
}
if (optConst) {
nll->constOptimizeTestStatistic(RooAbsArg::Activate,optConst>1) ;
}
if (doStripDisconnected) {
delete cPars ;
}
if (doOffset) {
nll->enableOffsetting(kTRUE) ;
}
return nll ;
}
RooFitResult* RooAbsPdf::fitTo(RooAbsData& data, const RooCmdArg& arg1, const RooCmdArg& arg2, const RooCmdArg& arg3, const RooCmdArg& arg4,
const RooCmdArg& arg5, const RooCmdArg& arg6, const RooCmdArg& arg7, const RooCmdArg& arg8)
{
RooLinkedList l ;
l.Add((TObject*)&arg1) ; l.Add((TObject*)&arg2) ;
l.Add((TObject*)&arg3) ; l.Add((TObject*)&arg4) ;
l.Add((TObject*)&arg5) ; l.Add((TObject*)&arg6) ;
l.Add((TObject*)&arg7) ; l.Add((TObject*)&arg8) ;
return fitTo(data,l) ;
}
RooFitResult* RooAbsPdf::fitTo(RooAbsData& data, const RooLinkedList& cmdList)
{
RooCmdConfig pc(Form("RooAbsPdf::fitTo(%s)",GetName())) ;
RooLinkedList fitCmdList(cmdList) ;
RooLinkedList nllCmdList = pc.filterCmdList(fitCmdList,"ProjectedObservables,Extended,Range,RangeWithName,SumCoefRange,NumCPU,SplitRange,Constrained,Constrain,ExternalConstraints,CloneData,GlobalObservables,GlobalObservablesTag,OffsetLikelihood") ;
pc.defineString("fitOpt","FitOptions",0,"") ;
pc.defineInt("optConst","Optimize",0,2) ;
pc.defineInt("verbose","Verbose",0,0) ;
pc.defineInt("doSave","Save",0,0) ;
pc.defineInt("doTimer","Timer",0,0) ;
pc.defineInt("plevel","PrintLevel",0,1) ;
pc.defineInt("strat","Strategy",0,1) ;
pc.defineInt("initHesse","InitialHesse",0,0) ;
pc.defineInt("hesse","Hesse",0,1) ;
pc.defineInt("minos","Minos",0,0) ;
pc.defineInt("ext","Extended",0,2) ;
pc.defineInt("numcpu","NumCPU",0,1) ;
pc.defineInt("numee","PrintEvalErrors",0,10) ;
pc.defineInt("doEEWall","EvalErrorWall",0,1) ;
pc.defineInt("doWarn","Warnings",0,1) ;
pc.defineInt("doSumW2","SumW2Error",0,-1) ;
pc.defineInt("doOffset","OffsetLikelihood",0,0) ;
pc.defineString("mintype","Minimizer",0,"OldMinuit") ;
pc.defineString("minalg","Minimizer",1,"minuit") ;
pc.defineObject("minosSet","Minos",0,0) ;
pc.defineSet("cPars","Constrain",0,0) ;
pc.defineSet("extCons","ExternalConstraints",0,0) ;
pc.defineMutex("FitOptions","Verbose") ;
pc.defineMutex("FitOptions","Save") ;
pc.defineMutex("FitOptions","Timer") ;
pc.defineMutex("FitOptions","Strategy") ;
pc.defineMutex("FitOptions","InitialHesse") ;
pc.defineMutex("FitOptions","Hesse") ;
pc.defineMutex("FitOptions","Minos") ;
pc.defineMutex("Range","RangeWithName") ;
pc.defineMutex("InitialHesse","Minimizer") ;
pc.process(fitCmdList) ;
if (!pc.ok(kTRUE)) {
return 0 ;
}
const char* fitOpt = pc.getString("fitOpt",0,kTRUE) ;
Int_t optConst = pc.getInt("optConst") ;
Int_t verbose = pc.getInt("verbose") ;
Int_t doSave = pc.getInt("doSave") ;
Int_t doTimer = pc.getInt("doTimer") ;
Int_t plevel = pc.getInt("plevel") ;
Int_t strat = pc.getInt("strat") ;
Int_t initHesse= pc.getInt("initHesse") ;
Int_t hesse = pc.getInt("hesse") ;
Int_t minos = pc.getInt("minos") ;
Int_t numee = pc.getInt("numee") ;
Int_t doEEWall = pc.getInt("doEEWall") ;
Int_t doWarn = pc.getInt("doWarn") ;
Int_t doSumW2 = pc.getInt("doSumW2") ;
const RooArgSet* minosSet = static_cast<RooArgSet*>(pc.getObject("minosSet")) ;
#ifdef __ROOFIT_NOROOMINIMIZER
const char* minType =0 ;
#else
const char* minType = pc.getString("mintype","OldMinuit") ;
const char* minAlg = pc.getString("minalg","minuit") ;
#endif
Bool_t weightedData = data.isNonPoissonWeighted() ;
if (weightedData && doSumW2==-1) {
coutW(InputArguments) << "RooAbsPdf::fitTo(" << GetName() << ") WARNING: a likelihood fit is request of what appears to be weighted data. " << endl
<< " While the estimated values of the parameters will always be calculated taking the weights into account, " << endl
<< " there are multiple ways to estimate the errors on these parameter values. You are advised to make an " << endl
<< " explicit choice on the error calculation: " << endl
<< " - Either provide SumW2Error(kTRUE), to calculate a sum-of-weights corrected HESSE error matrix " << endl
<< " (error will be proportional to the number of events)" << endl
<< " - Or provide SumW2Error(kFALSE), to return errors from original HESSE error matrix" << endl
<< " (which will be proportional to the sum of the weights)" << endl
<< " If you want the errors to reflect the information contained in the provided dataset, choose kTRUE. " << endl
<< " If you want the errors to reflect the precision you would be able to obtain with an unweighted dataset " << endl
<< " with 'sum-of-weights' events, choose kFALSE." << endl ;
}
if (doSumW2==1 && minos) {
coutW(InputArguments) << "RooAbsPdf::fitTo(" << GetName() << ") WARNING: sum-of-weights correction does not apply to MINOS errors" << endl ;
}
RooAbsReal* nll = createNLL(data,nllCmdList) ;
RooFitResult *ret = 0 ;
if (string(minType)!="OldMinuit") {
#ifndef __ROOFIT_NOROOMINIMIZER
RooMinimizer m(*nll) ;
m.setMinimizerType(minType) ;
m.setEvalErrorWall(doEEWall) ;
if (doWarn==0) {
}
m.setPrintEvalErrors(numee) ;
if (plevel!=1) {
m.setPrintLevel(plevel) ;
}
if (optConst) {
m.optimizeConst(optConst) ;
}
if (fitOpt) {
ret = m.fit(fitOpt) ;
} else {
if (verbose) {
m.setVerbose(1) ;
}
if (doTimer) {
m.setProfile(1) ;
}
if (strat!=1) {
m.setStrategy(strat) ;
}
if (initHesse) {
m.hesse() ;
}
m.minimize(minType,minAlg) ;
if (hesse) {
m.hesse() ;
}
if (doSumW2==1 && m.getNPar()>0) {
RooArgSet* comps = nll->getComponents();
vector<RooNLLVar*> nllComponents;
nllComponents.reserve(comps->getSize());
TIterator* citer = comps->createIterator();
RooAbsArg* arg;
while ((arg=(RooAbsArg*)citer->Next())) {
RooNLLVar* nllComp = dynamic_cast<RooNLLVar*>(arg);
if (!nllComp) continue;
nllComponents.push_back(nllComp);
}
delete citer;
delete comps;
RooFitResult* rw = m.save();
for (vector<RooNLLVar*>::iterator it = nllComponents.begin(); nllComponents.end() != it; ++it) {
(*it)->applyWeightSquared(kTRUE);
}
coutI(Fitting) << "RooAbsPdf::fitTo(" << GetName() << ") Calculating sum-of-weights-squared correction matrix for covariance matrix" << endl ;
m.hesse();
RooFitResult* rw2 = m.save();
for (vector<RooNLLVar*>::iterator it = nllComponents.begin(); nllComponents.end() != it; ++it) {
(*it)->applyWeightSquared(kFALSE);
}
const TMatrixDSym& matV = rw->covarianceMatrix();
TMatrixDSym matC = rw2->covarianceMatrix();
using ROOT::Math::CholeskyDecompGenDim;
CholeskyDecompGenDim<Double_t> decomp(matC.GetNrows(), matC);
if (!decomp) {
coutE(Fitting) << "RooAbsPdf::fitTo(" << GetName()
<< ") ERROR: Cannot apply sum-of-weights correction to covariance matrix: correction matrix calculated with weight-squared is singular" <<endl ;
} else {
decomp.Invert(matC);
for (int i = 0; i < matC.GetNrows(); ++i)
for (int j = 0; j < i; ++j) matC(j, i) = matC(i, j);
matC.Similarity(matV);
m.applyCovarianceMatrix(matC);
}
delete rw;
delete rw2;
}
if (minos) {
if (minosSet) {
m.minos(*minosSet) ;
} else {
m.minos() ;
}
}
if (doSave) {
string name = Form("fitresult_%s_%s",GetName(),data.GetName()) ;
string title = Form("Result of fit of p.d.f. %s to dataset %s",GetName(),data.GetName()) ;
ret = m.save(name.c_str(),title.c_str()) ;
}
}
if (optConst) {
m.optimizeConst(0) ;
}
#endif
} else {
RooMinuit m(*nll) ;
m.setEvalErrorWall(doEEWall) ;
if (doWarn==0) {
m.setNoWarn() ;
}
m.setPrintEvalErrors(numee) ;
if (plevel!=1) {
m.setPrintLevel(plevel) ;
}
if (optConst) {
m.optimizeConst(optConst) ;
}
if (fitOpt) {
ret = m.fit(fitOpt) ;
} else {
if (verbose) {
m.setVerbose(1) ;
}
if (doTimer) {
m.setProfile(1) ;
}
if (strat!=1) {
m.setStrategy(strat) ;
}
if (initHesse) {
m.hesse() ;
}
m.migrad() ;
if (hesse) {
m.hesse() ;
}
if (doSumW2==1 && m.getNPar()>0) {
list<RooNLLVar*> nllComponents ;
RooArgSet* comps = nll->getComponents() ;
RooAbsArg* arg ;
TIterator* citer = comps->createIterator() ;
while((arg=(RooAbsArg*)citer->Next())) {
RooNLLVar* nllComp = dynamic_cast<RooNLLVar*>(arg) ;
if (nllComp) {
nllComponents.push_back(nllComp) ;
}
}
delete citer ;
delete comps ;
RooFitResult* rw = m.save() ;
for (list<RooNLLVar*>::iterator iter1=nllComponents.begin() ; iter1!=nllComponents.end() ; iter1++) {
(*iter1)->applyWeightSquared(kTRUE) ;
}
coutI(Fitting) << "RooAbsPdf::fitTo(" << GetName() << ") Calculating sum-of-weights-squared correction matrix for covariance matrix" << endl ;
m.hesse() ;
RooFitResult* rw2 = m.save() ;
for (list<RooNLLVar*>::iterator iter2=nllComponents.begin() ; iter2!=nllComponents.end() ; iter2++) {
(*iter2)->applyWeightSquared(kFALSE) ;
}
const TMatrixDSym& matV = rw->covarianceMatrix();
TMatrixDSym matC = rw2->covarianceMatrix();
using ROOT::Math::CholeskyDecompGenDim;
CholeskyDecompGenDim<Double_t> decomp(matC.GetNrows(), matC);
if (!decomp) {
coutE(Fitting) << "RooAbsPdf::fitTo(" << GetName()
<< ") ERROR: Cannot apply sum-of-weights correction to covariance matrix: correction matrix calculated with weight-squared is singular" <<endl ;
} else {
decomp.Invert(matC);
for (int i = 0; i < matC.GetNrows(); ++i)
for (int j = 0; j < i; ++j) matC(j, i) = matC(i, j);
matC.Similarity(matV);
m.applyCovarianceMatrix(matC);
}
delete rw ;
delete rw2 ;
}
if (minos) {
if (minosSet) {
m.minos(*minosSet) ;
} else {
m.minos() ;
}
}
if (doSave) {
string name = Form("fitresult_%s_%s",GetName(),data.GetName()) ;
string title = Form("Result of fit of p.d.f. %s to dataset %s",GetName(),data.GetName()) ;
ret = m.save(name.c_str(),title.c_str()) ;
}
}
if (optConst) {
m.optimizeConst(0) ;
}
}
delete nll ;
return ret ;
}
RooFitResult* RooAbsPdf::chi2FitTo(RooDataHist& data, const RooLinkedList& cmdList)
{
RooCmdConfig pc(Form("RooAbsPdf::chi2FitTo(%s)",GetName())) ;
RooLinkedList fitCmdList(cmdList) ;
RooLinkedList chi2CmdList = pc.filterCmdList(fitCmdList,"Range,RangeWithName,NumCPU,Optimize,ProjectedObservables,AddCoefRange,SplitRange,DataError,Extended") ;
RooAbsReal* chi2 = createChi2(data,chi2CmdList) ;
RooFitResult* ret = chi2FitDriver(*chi2,fitCmdList) ;
delete chi2 ;
return ret ;
}
RooAbsReal* RooAbsPdf::createChi2(RooDataHist& data, const RooCmdArg& arg1, const RooCmdArg& arg2,
const RooCmdArg& arg3, const RooCmdArg& arg4, const RooCmdArg& arg5,
const RooCmdArg& arg6, const RooCmdArg& arg7, const RooCmdArg& arg8)
{
RooLinkedList cmdList ;
cmdList.Add((TObject*)&arg1) ; cmdList.Add((TObject*)&arg2) ;
cmdList.Add((TObject*)&arg3) ; cmdList.Add((TObject*)&arg4) ;
cmdList.Add((TObject*)&arg5) ; cmdList.Add((TObject*)&arg6) ;
cmdList.Add((TObject*)&arg7) ; cmdList.Add((TObject*)&arg8) ;
RooCmdConfig pc(Form("RooAbsPdf::createChi2(%s)",GetName())) ;
pc.defineString("rangeName","RangeWithName",0,"",kTRUE) ;
pc.allowUndefined(kTRUE) ;
pc.process(cmdList) ;
if (!pc.ok(kTRUE)) {
return 0 ;
}
const char* rangeName = pc.getString("rangeName",0,kTRUE) ;
RooAbsReal::setEvalErrorLoggingMode(RooAbsReal::CollectErrors) ;
RooAbsReal* chi2 ;
string baseName = Form("chi2_%s_%s",GetName(),data.GetName()) ;
if (!rangeName || strchr(rangeName,',')==0) {
chi2 = new RooChi2Var(baseName.c_str(),baseName.c_str(),*this,data,arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8) ;
} else {
const RooCmdArg* rarg(0) ;
string rcmd = "RangeWithName" ;
if (arg1.GetName()==rcmd) rarg = &arg1 ;
if (arg2.GetName()==rcmd) rarg = &arg2 ;
if (arg3.GetName()==rcmd) rarg = &arg3 ;
if (arg4.GetName()==rcmd) rarg = &arg4 ;
if (arg5.GetName()==rcmd) rarg = &arg5 ;
if (arg6.GetName()==rcmd) rarg = &arg6 ;
if (arg7.GetName()==rcmd) rarg = &arg7 ;
if (arg8.GetName()==rcmd) rarg = &arg8 ;
RooArgList chi2List ;
const size_t bufSize = strlen(rangeName)+1;
char* buf = new char[bufSize] ;
strlcpy(buf,rangeName,bufSize) ;
char* token = strtok(buf,",") ;
while(token) {
RooCmdArg subRangeCmd = RooFit::Range(token) ;
RooAbsReal* chi2Comp = new RooChi2Var(Form("%s_%s",baseName.c_str(),token),"chi^2",*this,data,
&arg1==rarg?subRangeCmd:arg1,&arg2==rarg?subRangeCmd:arg2,
&arg3==rarg?subRangeCmd:arg3,&arg4==rarg?subRangeCmd:arg4,
&arg5==rarg?subRangeCmd:arg5,&arg6==rarg?subRangeCmd:arg6,
&arg7==rarg?subRangeCmd:arg7,&arg8==rarg?subRangeCmd:arg8) ;
chi2List.add(*chi2Comp) ;
token = strtok(0,",") ;
}
delete[] buf ;
chi2 = new RooAddition(baseName.c_str(),"chi^2",chi2List,kTRUE) ;
}
RooAbsReal::setEvalErrorLoggingMode(RooAbsReal::PrintErrors) ;
return chi2 ;
}
RooAbsReal* RooAbsPdf::createChi2(RooDataSet& data, const RooLinkedList& cmdList)
{
RooCmdConfig pc(Form("RooAbsPdf::fitTo(%s)",GetName())) ;
pc.defineInt("integrate","Integrate",0,0) ;
pc.defineObject("yvar","YVar",0,0) ;
pc.process(cmdList) ;
if (!pc.ok(kTRUE)) {
return 0 ;
}
Bool_t integrate = pc.getInt("integrate") ;
RooRealVar* yvar = (RooRealVar*) pc.getObject("yvar") ;
string name = Form("chi2_%s_%s",GetName(),data.GetName()) ;
if (yvar) {
return new RooXYChi2Var(name.c_str(),name.c_str(),*this,data,*yvar,integrate) ;
} else {
return new RooXYChi2Var(name.c_str(),name.c_str(),*this,data,integrate) ;
}
}
void RooAbsPdf::printValue(ostream& os) const
{
getVal() ;
if (_norm) {
os << evaluate() << "/" << _norm->getVal() ;
} else {
os << evaluate() ;
}
}
void RooAbsPdf::printMultiline(ostream& os, Int_t contents, Bool_t verbose, TString indent) const
{
RooAbsReal::printMultiline(os,contents,verbose,indent);
os << indent << "--- RooAbsPdf ---" << endl;
os << indent << "Cached value = " << _value << endl ;
if (_norm) {
os << indent << " Normalization integral: " << endl ;
TString moreIndent(indent) ; moreIndent.Append(" ") ;
_norm->printStream(os,kName|kAddress|kTitle|kValue|kArgs,kSingleLine,moreIndent.Data()) ;
}
}
RooAbsGenContext* RooAbsPdf::binnedGenContext(const RooArgSet &vars, Bool_t verbose) const
{
return new RooBinnedGenContext(*this,vars,0,0,verbose) ;
}
RooAbsGenContext* RooAbsPdf::genContext(const RooArgSet &vars, const RooDataSet *prototype,
const RooArgSet* auxProto, Bool_t verbose) const
{
return new RooGenContext(*this,vars,prototype,auxProto,verbose) ;
}
RooAbsGenContext* RooAbsPdf::autoGenContext(const RooArgSet &vars, const RooDataSet* prototype, const RooArgSet* auxProto,
Bool_t verbose, Bool_t autoBinned, const char* binnedTag) const
{
if (prototype || (auxProto && auxProto->getSize()>0)) {
return genContext(vars,prototype,auxProto,verbose);
}
RooAbsGenContext *context(0) ;
if ( (autoBinned & isBinnedDistribution(vars)) || ( binnedTag && strlen(binnedTag) && (getAttribute(binnedTag)||string(binnedTag)=="*"))) {
context = binnedGenContext(vars,verbose) ;
} else {
context= genContext(vars,0,0,verbose);
}
return context ;
}
RooDataSet *RooAbsPdf::generate(const RooArgSet& whatVars, Int_t nEvents, const RooCmdArg& arg1,
const RooCmdArg& arg2, const RooCmdArg& arg3,const RooCmdArg& arg4, const RooCmdArg& arg5)
{
return generate(whatVars,RooFit::NumEvents(nEvents),arg1,arg2,arg3,arg4,arg5) ;
}
RooDataSet *RooAbsPdf::generate(const RooArgSet& whatVars, const RooCmdArg& arg1,const RooCmdArg& arg2,
const RooCmdArg& arg3,const RooCmdArg& arg4, const RooCmdArg& arg5,const RooCmdArg& arg6)
{
RooCmdConfig pc(Form("RooAbsPdf::generate(%s)",GetName())) ;
pc.defineObject("proto","PrototypeData",0,0) ;
pc.defineString("dsetName","Name",0,"") ;
pc.defineInt("randProto","PrototypeData",0,0) ;
pc.defineInt("resampleProto","PrototypeData",1,0) ;
pc.defineInt("verbose","Verbose",0,0) ;
pc.defineInt("extended","Extended",0,0) ;
pc.defineInt("nEvents","NumEvents",0,0) ;
pc.defineInt("autoBinned","AutoBinned",0,1) ;
pc.defineInt("expectedData","ExpectedData",0,0) ;
pc.defineDouble("nEventsD","NumEventsD",0,-1.) ;
pc.defineString("binnedTag","GenBinned",0,"") ;
pc.defineMutex("GenBinned","ProtoData") ;
pc.process(arg1,arg2,arg3,arg4,arg5,arg6) ;
if (!pc.ok(kTRUE)) {
return 0 ;
}
RooDataSet* protoData = static_cast<RooDataSet*>(pc.getObject("proto",0)) ;
const char* dsetName = pc.getString("dsetName") ;
Bool_t verbose = pc.getInt("verbose") ;
Bool_t randProto = pc.getInt("randProto") ;
Bool_t resampleProto = pc.getInt("resampleProto") ;
Bool_t extended = pc.getInt("extended") ;
Bool_t autoBinned = pc.getInt("autoBinned") ;
const char* binnedTag = pc.getString("binnedTag") ;
Int_t nEventsI = pc.getInt("nEvents") ;
Double_t nEventsD = pc.getInt("nEventsD") ;
Bool_t expectedData = pc.getInt("expectedData") ;
Double_t nEvents = (nEventsD>0) ? nEventsD : Double_t(nEventsI);
if (expectedData) {
binnedTag="*" ;
}
if (extended) {
if (nEvents == 0) nEvents = expectedEvents(&whatVars);
} else if (nEvents==0) {
cxcoutI(Generation) << "No number of events specified , number of events generated is "
<< GetName() << "::expectedEvents() = " << expectedEvents(&whatVars)<< endl ;
}
if (extended && protoData && !randProto) {
cxcoutI(Generation) << "WARNING Using generator option Extended() (Poisson distribution of #events) together "
<< "with a prototype dataset implies incomplete sampling or oversampling of proto data. "
<< "Set randomize flag in ProtoData() option to randomize prototype dataset order and thus "
<< "to randomize the set of over/undersampled prototype events for each generation cycle." << endl ;
}
RooDataSet* data ;
if (protoData) {
data = generate(whatVars,*protoData,Int_t(nEvents),verbose,randProto,resampleProto) ;
} else {
data = generate(whatVars,nEvents,verbose,autoBinned,binnedTag,expectedData, extended) ;
}
if (dsetName && strlen(dsetName)>0) {
data->SetName(dsetName) ;
}
return data ;
}
RooAbsPdf::GenSpec* RooAbsPdf::prepareMultiGen(const RooArgSet &whatVars,
const RooCmdArg& arg1,const RooCmdArg& arg2,
const RooCmdArg& arg3,const RooCmdArg& arg4,
const RooCmdArg& arg5,const RooCmdArg& arg6)
{
RooCmdConfig pc(Form("RooAbsPdf::generate(%s)",GetName())) ;
pc.defineObject("proto","PrototypeData",0,0) ;
pc.defineString("dsetName","Name",0,"") ;
pc.defineInt("randProto","PrototypeData",0,0) ;
pc.defineInt("resampleProto","PrototypeData",1,0) ;
pc.defineInt("verbose","Verbose",0,0) ;
pc.defineInt("extended","Extended",0,0) ;
pc.defineInt("nEvents","NumEvents",0,0) ;
pc.defineInt("autoBinned","AutoBinned",0,1) ;
pc.defineString("binnedTag","GenBinned",0,"") ;
pc.defineMutex("GenBinned","ProtoData") ;
pc.process(arg1,arg2,arg3,arg4,arg5,arg6) ;
if (!pc.ok(kTRUE)) {
return 0 ;
}
RooDataSet* protoData = static_cast<RooDataSet*>(pc.getObject("proto",0)) ;
const char* dsetName = pc.getString("dsetName") ;
Int_t nEvents = pc.getInt("nEvents") ;
Bool_t verbose = pc.getInt("verbose") ;
Bool_t randProto = pc.getInt("randProto") ;
Bool_t resampleProto = pc.getInt("resampleProto") ;
Bool_t extended = pc.getInt("extended") ;
Bool_t autoBinned = pc.getInt("autoBinned") ;
const char* binnedTag = pc.getString("binnedTag") ;
RooAbsGenContext* cx = autoGenContext(whatVars,protoData,0,verbose,autoBinned,binnedTag) ;
return new GenSpec(cx,whatVars,protoData,nEvents,extended,randProto,resampleProto,dsetName) ;
}
RooDataSet *RooAbsPdf::generate(RooAbsPdf::GenSpec& spec) const
{
Double_t nEvt = spec._nGen == 0 ? expectedEvents(spec._whatVars) : spec._nGen;
RooDataSet* ret = generate(*spec._genContext,spec._whatVars,spec._protoData, nEvt,kFALSE,spec._randProto,spec._resampleProto,
spec._init,spec._extended) ;
spec._init = kTRUE ;
return ret ;
}
RooDataSet *RooAbsPdf::generate(const RooArgSet &whatVars, Double_t nEvents, Bool_t verbose, Bool_t autoBinned, const char* binnedTag, Bool_t expectedData, Bool_t extended) const
{
if (nEvents==0 && extendMode()==CanNotBeExtended) {
return new RooDataSet("emptyData","emptyData",whatVars) ;
}
RooAbsGenContext *context = autoGenContext(whatVars,0,0,verbose,autoBinned,binnedTag) ;
if (expectedData) {
context->setExpectedData(kTRUE) ;
}
RooDataSet *generated = 0;
if(0 != context && context->isValid()) {
generated= context->generate(nEvents, kFALSE, extended);
}
else {
coutE(Generation) << "RooAbsPdf::generate(" << GetName() << ") cannot create a valid context" << endl;
}
if(0 != context) delete context;
return generated;
}
RooDataSet *RooAbsPdf::generate(RooAbsGenContext& context, const RooArgSet &whatVars, const RooDataSet *prototype,
Double_t nEvents, Bool_t , Bool_t randProtoOrder, Bool_t resampleProto,
Bool_t skipInit, Bool_t extended) const
{
if (nEvents==0 && (prototype==0 || prototype->numEntries()==0)) {
return new RooDataSet("emptyData","emptyData",whatVars) ;
}
RooDataSet *generated = 0;
if (resampleProto) {
randProtoOrder=kTRUE ;
}
if (randProtoOrder && prototype && prototype->numEntries()!=nEvents) {
coutI(Generation) << "RooAbsPdf::generate (Re)randomizing event order in prototype dataset (Nevt=" << nEvents << ")" << endl ;
Int_t* newOrder = randomizeProtoOrder(prototype->numEntries(),Int_t(nEvents),resampleProto) ;
context.setProtoDataOrder(newOrder) ;
delete[] newOrder ;
}
if(context.isValid()) {
generated= context.generate(nEvents,skipInit,extended);
}
else {
coutE(Generation) << "RooAbsPdf::generate(" << GetName() << ") do not have a valid generator context" << endl;
}
return generated;
}
RooDataSet *RooAbsPdf::generate(const RooArgSet &whatVars, const RooDataSet& prototype,
Int_t nEvents, Bool_t verbose, Bool_t randProtoOrder, Bool_t resampleProto) const
{
RooAbsGenContext *context= genContext(whatVars,&prototype,0,verbose);
if (context) {
RooDataSet* data = generate(*context,whatVars,&prototype,nEvents,verbose,randProtoOrder,resampleProto) ;
delete context ;
return data ;
} else {
coutE(Generation) << "RooAbsPdf::generate(" << GetName() << ") ERROR creating generator context" << endl ;
return 0 ;
}
}
Int_t* RooAbsPdf::randomizeProtoOrder(Int_t nProto, Int_t, Bool_t resampleProto) const
{
RooLinkedList l ;
Int_t i ;
for (i=0 ; i<nProto ; i++) {
l.Add(new RooInt(i)) ;
}
Int_t* lut = new Int_t[nProto] ;
if (!resampleProto) {
for (i=0 ; i<nProto ; i++) {
Int_t iran = RooRandom::integer(nProto-i) ;
RooInt* sample = (RooInt*) l.At(iran) ;
lut[i] = *sample ;
l.Remove(sample) ;
delete sample ;
}
} else {
for (i=0 ; i<nProto ; i++) {
lut[i] = RooRandom::integer(nProto);
}
}
return lut ;
}
Int_t RooAbsPdf::getGenerator(const RooArgSet &, RooArgSet &, Bool_t ) const
{
return 0 ;
}
void RooAbsPdf::initGenerator(Int_t )
{
}
void RooAbsPdf::generateEvent(Int_t )
{
}
Bool_t RooAbsPdf::isDirectGenSafe(const RooAbsArg& arg) const
{
if (!findServer(arg.GetName())) return kFALSE ;
TIterator* sIter = serverIterator() ;
const RooAbsArg *server = 0;
while((server=(const RooAbsArg*)sIter->Next())) {
if(server == &arg) continue;
if(server->dependsOn(arg)) {
delete sIter ;
return kFALSE ;
}
}
delete sIter ;
return kTRUE ;
}
RooDataHist *RooAbsPdf::generateBinned(const RooArgSet& whatVars, Double_t nEvents, const RooCmdArg& arg1,
const RooCmdArg& arg2, const RooCmdArg& arg3,const RooCmdArg& arg4, const RooCmdArg& arg5)
{
return generateBinned(whatVars,RooFit::NumEvents(nEvents),arg1,arg2,arg3,arg4,arg5) ;
}
RooDataHist *RooAbsPdf::generateBinned(const RooArgSet& whatVars, const RooCmdArg& arg1,const RooCmdArg& arg2,
const RooCmdArg& arg3,const RooCmdArg& arg4, const RooCmdArg& arg5,const RooCmdArg& arg6)
{
RooCmdConfig pc(Form("RooAbsPdf::generate(%s)",GetName())) ;
pc.defineString("dsetName","Name",0,"") ;
pc.defineInt("verbose","Verbose",0,0) ;
pc.defineInt("extended","Extended",0,0) ;
pc.defineInt("nEvents","NumEvents",0,0) ;
pc.defineDouble("nEventsD","NumEventsD",0,-1.) ;
pc.defineInt("expectedData","ExpectedData",0,0) ;
pc.process(arg1,arg2,arg3,arg4,arg5,arg6) ;
if (!pc.ok(kTRUE)) {
return 0 ;
}
Double_t nEvents = pc.getDouble("nEventsD") ;
if (nEvents<0) {
nEvents = pc.getInt("nEvents") ;
}
Bool_t extended = pc.getInt("extended") ;
Bool_t expectedData = pc.getInt("expectedData") ;
const char* dsetName = pc.getString("dsetName") ;
if (extended) {
nEvents = (nEvents==0 ? expectedEvents(&whatVars) :nEvents) ;
cxcoutI(Generation) << " Extended mode active, number of events generated (" << nEvents << ") is Poisson fluctuation on "
<< GetName() << "::expectedEvents() = " << nEvents << endl ;
if (nEvents==0) {
return 0 ;
}
} else if (nEvents==0) {
cxcoutI(Generation) << "No number of events specified , number of events generated is "
<< GetName() << "::expectedEvents() = " << expectedEvents(&whatVars)<< endl ;
}
RooDataHist* data = generateBinned(whatVars,nEvents,expectedData,extended) ;
if (dsetName && strlen(dsetName)>0) {
data->SetName(dsetName) ;
}
return data ;
}
RooDataHist *RooAbsPdf::generateBinned(const RooArgSet &whatVars, Double_t nEvents, Bool_t expectedData, Bool_t extended) const
{
RooDataHist* hist = new RooDataHist("genData","genData",whatVars) ;
if (nEvents<=0) {
if (!canBeExtended()) {
coutE(InputArguments) << "RooAbsPdf::generateBinned(" << GetName() << ") ERROR: No event count provided and p.d.f does not provide expected number of events" << endl ;
delete hist ;
return 0 ;
} else {
if (expectedData || extended) {
nEvents = expectedEvents(&whatVars) ;
} else {
nEvents = Int_t(expectedEvents(&whatVars)+0.5) ;
}
}
}
fillDataHist(hist,&whatVars,1,kTRUE) ;
vector<int> histOut(hist->numEntries()) ;
Double_t histMax(-1) ;
Int_t histOutSum(0) ;
for (int i=0 ; i<hist->numEntries() ; i++) {
hist->get(i) ;
if (expectedData) {
Double_t w=hist->weight()*nEvents ;
hist->set(w,sqrt(w)) ;
} else if (extended) {
Double_t w = RooRandom::randomGenerator()->Poisson(hist->weight()*nEvents) ;
hist->set(w,sqrt(w)) ;
} else {
if (hist->weight()>histMax) {
histMax = hist->weight() ;
}
histOut[i] = RooRandom::randomGenerator()->Poisson(hist->weight()*nEvents) ;
histOutSum += histOut[i] ;
}
}
if (!expectedData && !extended) {
Int_t nEvtExtra = abs(Int_t(nEvents)-histOutSum) ;
Int_t wgt = (histOutSum>nEvents) ? -1 : 1 ;
while(nEvtExtra>0) {
Int_t ibinRand = RooRandom::randomGenerator()->Integer(hist->numEntries()) ;
hist->get(ibinRand) ;
Double_t ranY = RooRandom::randomGenerator()->Uniform(histMax) ;
if (ranY<hist->weight()) {
if (wgt==1) {
histOut[ibinRand]++ ;
} else {
if (histOut[ibinRand]>0) {
histOut[ibinRand]-- ;
} else {
continue ;
}
}
nEvtExtra-- ;
}
}
for (int i=0 ; i<hist->numEntries() ; i++) {
hist->get(i) ;
hist->set(histOut[i],sqrt(1.0*histOut[i])) ;
}
} else if (expectedData) {
Double_t corr = nEvents/hist->sumEntries() ;
for (int i=0 ; i<hist->numEntries() ; i++) {
hist->get(i) ;
hist->set(hist->weight()*corr,sqrt(hist->weight()*corr)) ;
}
}
return hist;
}
RooDataSet* RooAbsPdf::generateSimGlobal(const RooArgSet& whatVars, Int_t nEvents)
{
return generate(whatVars,nEvents) ;
}
RooPlot* RooAbsPdf::plotOn(RooPlot* frame, RooLinkedList& cmdList) const
{
RooCmdArg* plotRange(0) ;
RooCmdArg* normRange2(0) ;
if (getStringAttribute("fitrange") && !cmdList.FindObject("Range") &&
!cmdList.FindObject("RangeWithName")) {
plotRange = (RooCmdArg*) RooFit::Range(getStringAttribute("fitrange")).Clone() ;
cmdList.Add(plotRange) ;
}
if (getStringAttribute("fitrange") && !cmdList.FindObject("NormRange")) {
normRange2 = (RooCmdArg*) RooFit::NormRange(getStringAttribute("fitrange")).Clone() ;
cmdList.Add(normRange2) ;
}
if (plotRange || normRange2) {
coutI(Plotting) << "RooAbsPdf::plotOn(" << GetName() << ") p.d.f was fitted in range and no explicit "
<< (plotRange?"plot":"") << ((plotRange&&normRange2)?",":"")
<< (normRange2?"norm":"") << " range was specified, using fit range as default" << endl ;
}
if (plotSanityChecks(frame)) return frame ;
RooCmdConfig pc(Form("RooAbsPdf::plotOn(%s)",GetName())) ;
pc.defineDouble("scaleFactor","Normalization",0,1.0) ;
pc.defineInt("scaleType","Normalization",0,Relative) ;
pc.defineObject("compSet","SelectCompSet",0) ;
pc.defineString("compSpec","SelectCompSpec",0) ;
pc.defineObject("asymCat","Asymmetry",0) ;
pc.defineDouble("rangeLo","Range",0,-999.) ;
pc.defineDouble("rangeHi","Range",1,-999.) ;
pc.defineString("rangeName","RangeWithName",0,"") ;
pc.defineString("normRangeName","NormRange",0,"") ;
pc.defineInt("rangeAdjustNorm","Range",0,0) ;
pc.defineInt("rangeWNAdjustNorm","RangeWithName",0,0) ;
pc.defineMutex("SelectCompSet","SelectCompSpec") ;
pc.defineMutex("Range","RangeWithName") ;
pc.allowUndefined() ;
pc.process(cmdList) ;
if (!pc.ok(kTRUE)) {
return frame ;
}
ScaleType stype = (ScaleType) pc.getInt("scaleType") ;
Double_t scaleFactor = pc.getDouble("scaleFactor") ;
const RooAbsCategoryLValue* asymCat = (const RooAbsCategoryLValue*) pc.getObject("asymCat") ;
const char* compSpec = pc.getString("compSpec") ;
const RooArgSet* compSet = (const RooArgSet*) pc.getObject("compSet") ;
Bool_t haveCompSel = ( (compSpec && strlen(compSpec)>0) || compSet) ;
TString nameSuffix ;
if (compSpec && strlen(compSpec)>0) {
nameSuffix.Append("_Comp[") ;
nameSuffix.Append(compSpec) ;
nameSuffix.Append("]") ;
} else if (compSet) {
nameSuffix.Append("_Comp[") ;
nameSuffix.Append(compSet->contentsString().c_str()) ;
nameSuffix.Append("]") ;
}
pc.stripCmdList(cmdList,"SelectCompSet,SelectCompSpec") ;
if (asymCat) {
RooCmdArg cnsuffix("CurveNameSuffix",0,0,0,0,nameSuffix.Data(),0,0,0) ;
cmdList.Add(&cnsuffix);
return RooAbsReal::plotOn(frame,cmdList) ;
}
Double_t nExpected(1) ;
if (stype==RelativeExpected) {
if (!canBeExtended()) {
coutE(Plotting) << "RooAbsPdf::plotOn(" << GetName()
<< "): ERROR the 'Expected' scale option can only be used on extendable PDFs" << endl ;
return frame ;
}
nExpected = expectedEvents(frame->getNormVars()) ;
}
if (stype != Raw) {
if (frame->getFitRangeNEvt() && stype==Relative) {
Bool_t hasCustomRange(kFALSE), adjustNorm(kFALSE) ;
list<pair<Double_t,Double_t> > rangeLim ;
if (pc.hasProcessed("Range")) {
Double_t rangeLo = pc.getDouble("rangeLo") ;
Double_t rangeHi = pc.getDouble("rangeHi") ;
rangeLim.push_back(make_pair(rangeLo,rangeHi)) ;
adjustNorm = pc.getInt("rangeAdjustNorm") ;
hasCustomRange = kTRUE ;
coutI(Plotting) << "RooAbsPdf::plotOn(" << GetName() << ") only plotting range ["
<< rangeLo << "," << rangeHi << "]" ;
if (!pc.hasProcessed("NormRange")) {
ccoutI(Plotting) << ", curve is normalized to data in " << (adjustNorm?"given":"full") << " given range" << endl ;
} else {
ccoutI(Plotting) << endl ;
}
nameSuffix.Append(Form("_Range[%f_%f]",rangeLo,rangeHi)) ;
} else if (pc.hasProcessed("RangeWithName")) {
char tmp[1024] ;
strlcpy(tmp,pc.getString("rangeName",0,kTRUE),1024) ;
char* rangeNameToken = strtok(tmp,",") ;
while(rangeNameToken) {
Double_t rangeLo = frame->getPlotVar()->getMin(rangeNameToken) ;
Double_t rangeHi = frame->getPlotVar()->getMax(rangeNameToken) ;
rangeLim.push_back(make_pair(rangeLo,rangeHi)) ;
rangeNameToken = strtok(0,",") ;
}
adjustNorm = pc.getInt("rangeWNAdjustNorm") ;
hasCustomRange = kTRUE ;
coutI(Plotting) << "RooAbsPdf::plotOn(" << GetName() << ") only plotting range '" << pc.getString("rangeName",0,kTRUE) << "'" ;
if (!pc.hasProcessed("NormRange")) {
ccoutI(Plotting) << ", curve is normalized to data in " << (adjustNorm?"given":"full") << " given range" << endl ;
} else {
ccoutI(Plotting) << endl ;
}
nameSuffix.Append(Form("_Range[%s]",pc.getString("rangeName"))) ;
}
if (pc.hasProcessed("NormRange")) {
char tmp[1024] ;
strlcpy(tmp,pc.getString("normRangeName",0,kTRUE),1024) ;
char* rangeNameToken = strtok(tmp,",") ;
rangeLim.clear() ;
while(rangeNameToken) {
Double_t rangeLo = frame->getPlotVar()->getMin(rangeNameToken) ;
Double_t rangeHi = frame->getPlotVar()->getMax(rangeNameToken) ;
rangeLim.push_back(make_pair(rangeLo,rangeHi)) ;
rangeNameToken = strtok(0,",") ;
}
adjustNorm = kTRUE ;
hasCustomRange = kTRUE ;
coutI(Plotting) << "RooAbsPdf::plotOn(" << GetName() << ") p.d.f. curve is normalized using explicit choice of ranges '" << pc.getString("normRangeName",0,kTRUE) << "'" << endl ;
nameSuffix.Append(Form("_NormRange[%s]",pc.getString("rangeName"))) ;
}
if (hasCustomRange && adjustNorm) {
Double_t rangeNevt(0) ;
list<pair<Double_t,Double_t> >::iterator riter = rangeLim.begin() ;
for (;riter!=rangeLim.end() ; ++riter) {
Double_t nevt= frame->getFitRangeNEvt(riter->first,riter->second) ;
rangeNevt += nevt ;
}
scaleFactor *= rangeNevt/nExpected ;
} else {
scaleFactor *= frame->getFitRangeNEvt()/nExpected ;
}
} else if (stype==RelativeExpected) {
scaleFactor *= nExpected ;
} else if (stype==NumEvent) {
scaleFactor /= nExpected ;
}
scaleFactor *= frame->getFitRangeBinW() ;
}
frame->updateNormVars(*frame->getPlotVar()) ;
RooCmdArg tmp = RooFit::Normalization(scaleFactor,Raw) ;
tmp.setInt(1,1) ;
cmdList.Add(&tmp) ;
if (haveCompSel) {
RooArgSet branchNodeSet ;
branchNodeServerList(&branchNodeSet) ;
TIterator* iter = branchNodeSet.createIterator() ;
RooAbsArg* arg ;
while((arg=(RooAbsArg*)iter->Next())) {
if (!dynamic_cast<RooAbsReal*>(arg)) {
branchNodeSet.remove(*arg) ;
}
}
delete iter ;
RooArgSet* dirSelNodes ;
if (compSet) {
dirSelNodes = (RooArgSet*) branchNodeSet.selectCommon(*compSet) ;
} else {
dirSelNodes = (RooArgSet*) branchNodeSet.selectByName(compSpec) ;
}
if (dirSelNodes->getSize()>0) {
coutI(Plotting) << "RooAbsPdf::plotOn(" << GetName() << ") directly selected PDF components: " << *dirSelNodes << endl ;
plotOnCompSelect(dirSelNodes) ;
} else {
if (compSet) {
coutE(Plotting) << "RooAbsPdf::plotOn(" << GetName() << ") ERROR: component selection set " << *compSet << " does not match any components of p.d.f." << endl ;
} else {
coutE(Plotting) << "RooAbsPdf::plotOn(" << GetName() << ") ERROR: component selection expression '" << compSpec << "' does not select any components of p.d.f." << endl ;
}
return 0 ;
}
delete dirSelNodes ;
}
RooCmdArg cnsuffix("CurveNameSuffix",0,0,0,0,nameSuffix.Data(),0,0,0) ;
cmdList.Add(&cnsuffix);
RooPlot* ret = RooAbsReal::plotOn(frame,cmdList) ;
if (haveCompSel) plotOnCompSelect(0) ;
if (plotRange) {
delete plotRange ;
}
if (normRange2) {
delete normRange2 ;
}
return ret ;
}
RooPlot* RooAbsPdf::plotOn(RooPlot *frame, PlotOpt o) const
{
if (plotSanityChecks(frame)) return frame ;
Double_t nExpected(1) ;
if (o.stype==RelativeExpected) {
if (!canBeExtended()) {
coutE(Plotting) << "RooAbsPdf::plotOn(" << GetName()
<< "): ERROR the 'Expected' scale option can only be used on extendable PDFs" << endl ;
return frame ;
}
nExpected = expectedEvents(frame->getNormVars()) ;
}
if (o.stype != Raw) {
if (frame->getFitRangeNEvt() && o.stype==Relative) {
o.scaleFactor *= frame->getFitRangeNEvt()/nExpected ;
} else if (o.stype==RelativeExpected) {
o.scaleFactor *= nExpected ;
} else if (o.stype==NumEvent) {
o.scaleFactor /= nExpected ;
}
o.scaleFactor *= frame->getFitRangeBinW() ;
}
frame->updateNormVars(*frame->getPlotVar()) ;
return RooAbsReal::plotOn(frame,o) ;
}
RooPlot* RooAbsPdf::paramOn(RooPlot* frame, const RooCmdArg& arg1, const RooCmdArg& arg2,
const RooCmdArg& arg3, const RooCmdArg& arg4, const RooCmdArg& arg5,
const RooCmdArg& arg6, const RooCmdArg& arg7, const RooCmdArg& arg8)
{
RooLinkedList cmdList;
cmdList.Add(const_cast<RooCmdArg*>(&arg1)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg2)) ;
cmdList.Add(const_cast<RooCmdArg*>(&arg3)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg4)) ;
cmdList.Add(const_cast<RooCmdArg*>(&arg5)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg6)) ;
cmdList.Add(const_cast<RooCmdArg*>(&arg7)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg8)) ;
RooCmdConfig pc(Form("RooAbsPdf::paramOn(%s)",GetName())) ;
pc.defineString("label","Label",0,"") ;
pc.defineDouble("xmin","Layout",0,0.50) ;
pc.defineDouble("xmax","Layout",1,0.99) ;
pc.defineInt("ymaxi","Layout",0,Int_t(0.95*10000)) ;
pc.defineInt("showc","ShowConstants",0,0) ;
pc.defineObject("params","Parameters",0,0) ;
pc.defineString("formatStr","Format",0,"NELU") ;
pc.defineInt("sigDigit","Format",0,2) ;
pc.defineInt("dummy","FormatArgs",0,0) ;
pc.defineMutex("Format","FormatArgs") ;
pc.process(cmdList) ;
if (!pc.ok(kTRUE)) {
return frame ;
}
const char* label = pc.getString("label") ;
Double_t xmin = pc.getDouble("xmin") ;
Double_t xmax = pc.getDouble("xmax") ;
Double_t ymax = pc.getInt("ymaxi") / 10000. ;
Int_t showc = pc.getInt("showc") ;
const char* formatStr = pc.getString("formatStr") ;
Int_t sigDigit = pc.getInt("sigDigit") ;
RooArgSet* params = static_cast<RooArgSet*>(pc.getObject("params")) ;
if (!params) {
params = getParameters(frame->getNormVars()) ;
if (pc.hasProcessed("FormatArgs")) {
const RooCmdArg* formatCmd = static_cast<RooCmdArg*>(cmdList.FindObject("FormatArgs")) ;
paramOn(frame,*params,showc,label,0,0,xmin,xmax,ymax,formatCmd) ;
} else {
paramOn(frame,*params,showc,label,sigDigit,formatStr,xmin,xmax,ymax) ;
}
delete params ;
} else {
RooArgSet* pdfParams = getParameters(frame->getNormVars()) ;
RooArgSet* selParams = static_cast<RooArgSet*>(pdfParams->selectCommon(*params)) ;
if (pc.hasProcessed("FormatArgs")) {
const RooCmdArg* formatCmd = static_cast<RooCmdArg*>(cmdList.FindObject("FormatArgs")) ;
paramOn(frame,*selParams,showc,label,0,0,xmin,xmax,ymax,formatCmd) ;
} else {
paramOn(frame,*selParams,showc,label,sigDigit,formatStr,xmin,xmax,ymax) ;
}
delete selParams ;
delete pdfParams ;
}
return frame ;
}
RooPlot* RooAbsPdf::paramOn(RooPlot* frame, const RooAbsData* data, const char *label,
Int_t sigDigits, Option_t *options, Double_t xmin,
Double_t xmax ,Double_t ymax)
{
RooArgSet* params = getParameters(data) ;
TString opts(options) ;
paramOn(frame,*params,opts.Contains("c"),label,sigDigits,options,xmin,xmax,ymax) ;
delete params ;
return frame ;
}
RooPlot* RooAbsPdf::paramOn(RooPlot* frame, const RooArgSet& params, Bool_t showConstants, const char *label,
Int_t sigDigits, Option_t *options, Double_t xmin,
Double_t xmax ,Double_t ymax, const RooCmdArg* formatCmd)
{
TString opts = options;
opts.ToLower();
Bool_t showLabel= (label != 0 && strlen(label) > 0);
TIterator* pIter = params.createIterator() ;
Double_t ymin(ymax), dy(0.06);
RooRealVar *var = 0;
while((var=(RooRealVar*)pIter->Next())) {
if(showConstants || !var->isConstant()) ymin-= dy;
}
if(showLabel) ymin-= dy;
TPaveText *box= new TPaveText(xmin,ymax,xmax,ymin,"BRNDC");
if(!box) return 0;
box->SetName(Form("%s_paramBox",GetName())) ;
box->SetFillColor(0);
box->SetBorderSize(1);
box->SetTextAlign(12);
box->SetTextSize(0.04F);
box->SetFillStyle(1001);
box->SetFillColor(0);
pIter->Reset() ;
while((var=(RooRealVar*)pIter->Next())) {
if(var->isConstant() && !showConstants) continue;
TString *formatted= options ? var->format(sigDigits, options) : var->format(*formatCmd) ;
box->AddText(formatted->Data());
delete formatted;
}
if(showLabel) box->AddText(label);
frame->addObject(box) ;
delete pIter ;
return frame ;
}
Double_t RooAbsPdf::expectedEvents(const RooArgSet*) const
{
return 0 ;
}
void RooAbsPdf::verboseEval(Int_t stat)
{
_verboseEval = stat ;
}
Int_t RooAbsPdf::verboseEval()
{
return _verboseEval ;
}
void RooAbsPdf::CacheElem::operModeHook(RooAbsArg::OperMode)
{
}
RooAbsPdf::CacheElem::~CacheElem()
{
if (_owner) {
RooAbsPdf* pdfOwner = static_cast<RooAbsPdf*>(_owner) ;
if (pdfOwner->_norm == _norm) {
pdfOwner->_norm = 0 ;
}
}
delete _norm ;
}
RooAbsPdf* RooAbsPdf::createProjection(const RooArgSet& iset)
{
TString name(GetName()) ;
name.Append("_Proj[") ;
if (iset.getSize()>0) {
TIterator* iter = iset.createIterator() ;
RooAbsArg* arg ;
Bool_t first(kTRUE) ;
while((arg=(RooAbsArg*)iter->Next())) {
if (first) {
first=kFALSE ;
} else {
name.Append(",") ;
}
name.Append(arg->GetName()) ;
}
delete iter ;
}
name.Append("]") ;
return new RooProjectedPdf(name.Data(),name.Data(),*this,iset) ;
}
RooAbsReal* RooAbsPdf::createCdf(const RooArgSet& iset, const RooArgSet& nset)
{
return createCdf(iset,RooFit::SupNormSet(nset)) ;
}
RooAbsReal* RooAbsPdf::createCdf(const RooArgSet& iset, const RooCmdArg& arg1, const RooCmdArg& arg2,
const RooCmdArg& arg3, const RooCmdArg& arg4, const RooCmdArg& arg5,
const RooCmdArg& arg6, const RooCmdArg& arg7, const RooCmdArg& arg8)
{
RooCmdConfig pc(Form("RooAbsReal::createCdf(%s)",GetName())) ;
pc.defineObject("supNormSet","SupNormSet",0,0) ;
pc.defineInt("numScanBins","ScanParameters",0,1000) ;
pc.defineInt("intOrder","ScanParameters",1,2) ;
pc.defineInt("doScanNum","ScanNumCdf",0,1) ;
pc.defineInt("doScanAll","ScanAllCdf",0,0) ;
pc.defineInt("doScanNon","ScanNoCdf",0,0) ;
pc.defineMutex("ScanNumCdf","ScanAllCdf","ScanNoCdf") ;
pc.process(arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8) ;
if (!pc.ok(kTRUE)) {
return 0 ;
}
const RooArgSet* snset = static_cast<const RooArgSet*>(pc.getObject("supNormSet",0)) ;
RooArgSet nset ;
if (snset) {
nset.add(*snset) ;
}
Int_t numScanBins = pc.getInt("numScanBins") ;
Int_t intOrder = pc.getInt("intOrder") ;
Int_t doScanNum = pc.getInt("doScanNum") ;
Int_t doScanAll = pc.getInt("doScanAll") ;
Int_t doScanNon = pc.getInt("doScanNon") ;
if (doScanNon) {
return createIntRI(iset,nset) ;
}
if (doScanAll) {
return createScanCdf(iset,nset,numScanBins,intOrder) ;
}
if (doScanNum) {
RooRealIntegral* tmp = (RooRealIntegral*) createIntegral(iset) ;
Int_t isNum= (tmp->numIntRealVars().getSize()>0) ;
delete tmp ;
if (isNum) {
coutI(NumIntegration) << "RooAbsPdf::createCdf(" << GetName() << ") integration over observable(s) " << iset << " involves numeric integration," << endl
<< " constructing cdf though numeric integration of sampled pdf in " << numScanBins << " bins and applying order "
<< intOrder << " interpolation on integrated histogram." << endl
<< " To override this choice of technique use argument ScanNone(), to change scan parameters use ScanParameters(nbins,order) argument" << endl ;
}
return isNum ? createScanCdf(iset,nset,numScanBins,intOrder) : createIntRI(iset,nset) ;
}
return 0 ;
}
RooAbsReal* RooAbsPdf::createScanCdf(const RooArgSet& iset, const RooArgSet& nset, Int_t numScanBins, Int_t intOrder)
{
string name = string(GetName()) + "_NUMCDF_" + integralNameSuffix(iset,&nset).Data() ;
RooRealVar* ivar = (RooRealVar*) iset.first() ;
ivar->setBins(numScanBins,"numcdf") ;
RooNumCdf* ret = new RooNumCdf(name.c_str(),name.c_str(),*this,*ivar,"numcdf") ;
ret->setInterpolationOrder(intOrder) ;
return ret ;
}
RooArgSet* RooAbsPdf::getAllConstraints(const RooArgSet& observables, RooArgSet& constrainedParams, Bool_t stripDisconnected) const
{
RooArgSet* ret = new RooArgSet("AllConstraints") ;
RooArgSet* comps = getComponents() ;
TIterator* iter = comps->createIterator() ;
RooAbsArg *arg ;
while((arg=(RooAbsArg*)iter->Next())) {
RooAbsPdf* pdf = dynamic_cast<RooAbsPdf*>(arg) ;
if (pdf && !ret->find(pdf->GetName())) {
RooArgSet* compRet = pdf->getConstraints(observables,constrainedParams,stripDisconnected) ;
if (compRet) {
ret->add(*compRet,kFALSE) ;
delete compRet ;
}
}
}
delete iter ;
delete comps ;
return ret ;
}
void RooAbsPdf::clearEvalError()
{
_evalError = kFALSE ;
}
Bool_t RooAbsPdf::evalError()
{
return _evalError ;
}
void RooAbsPdf::raiseEvalError()
{
_evalError = kTRUE ;
}
RooNumGenConfig* RooAbsPdf::defaultGeneratorConfig()
{
return &RooNumGenConfig::defaultConfig() ;
}
RooNumGenConfig* RooAbsPdf::specialGeneratorConfig() const
{
return _specGeneratorConfig ;
}
RooNumGenConfig* RooAbsPdf::specialGeneratorConfig(Bool_t createOnTheFly)
{
if (!_specGeneratorConfig && createOnTheFly) {
_specGeneratorConfig = new RooNumGenConfig(*defaultGeneratorConfig()) ;
}
return _specGeneratorConfig ;
}
const RooNumGenConfig* RooAbsPdf::getGeneratorConfig() const
{
const RooNumGenConfig* config = specialGeneratorConfig() ;
if (config) return config ;
return defaultGeneratorConfig() ;
}
void RooAbsPdf::setGeneratorConfig(const RooNumGenConfig& config)
{
if (_specGeneratorConfig) {
delete _specGeneratorConfig ;
}
_specGeneratorConfig = new RooNumGenConfig(config) ;
}
void RooAbsPdf::setGeneratorConfig()
{
if (_specGeneratorConfig) {
delete _specGeneratorConfig ;
}
_specGeneratorConfig = 0 ;
}
RooAbsPdf::GenSpec::~GenSpec()
{
delete _genContext ;
}
RooAbsPdf::GenSpec::GenSpec(RooAbsGenContext* context, const RooArgSet& whatVars, RooDataSet* protoData, Int_t nGen,
Bool_t extended, Bool_t randProto, Bool_t resampleProto, TString dsetName, Bool_t init) :
_genContext(context), _whatVars(whatVars), _protoData(protoData), _nGen(nGen), _extended(extended),
_randProto(randProto), _resampleProto(resampleProto), _dsetName(dsetName), _init(init)
{
}
void RooAbsPdf::setNormRange(const char* rangeName)
{
if (rangeName) {
_normRange = rangeName ;
} else {
_normRange.Clear() ;
}
if (_norm) {
_normMgr.sterilize() ;
_norm = 0 ;
}
}
void RooAbsPdf::setNormRangeOverride(const char* rangeName)
{
if (rangeName) {
_normRangeOverride = rangeName ;
} else {
_normRangeOverride.Clear() ;
}
if (_norm) {
_normMgr.sterilize() ;
_norm = 0 ;
}
}