#include <iostream>
#include <fstream>
#include "RooStats/HLFactory.h"
#include "TFile.h"
#include "TObject.h"
#include "TObjArray.h"
#include "TObjString.h"
#include "RooSimultaneous.h"
using namespace std;
ClassImp(RooStats::HLFactory) ;
using namespace RooStats;
using namespace RooFit;
HLFactory::HLFactory(const char *name,
const char *fileName,
bool isVerbose):
TNamed(name,name),
fComboCat(0),
fComboBkgPdf(0),
fComboSigBkgPdf(0),
fComboDataset(0),
fCombinationDone(false),
fVerbose(isVerbose),
fInclusionLevel(0),
fOwnWs(true){
TString wsName(name);
wsName+="_ws";
fWs = new RooWorkspace(wsName,true);
fSigBkgPdfNames.SetOwner();
fBkgPdfNames.SetOwner();
fDatasetsNames.SetOwner();
fReadFile(fileName);
}
HLFactory::HLFactory(const char* name,
RooWorkspace* externalWs,
bool isVerbose):
TNamed(name,name),
fComboCat(0),
fComboBkgPdf(0),
fComboSigBkgPdf(0),
fComboDataset(0),
fCombinationDone(false),
fVerbose(isVerbose),
fInclusionLevel(0),
fOwnWs(false){
fWs=externalWs;
fSigBkgPdfNames.SetOwner();
fBkgPdfNames.SetOwner();
fDatasetsNames.SetOwner();
}
HLFactory::HLFactory():
TNamed("hlfactory","hlfactory"),
fComboCat(0),
fComboBkgPdf(0),
fComboSigBkgPdf(0),
fComboDataset(0),
fCombinationDone(false),
fVerbose(false),
fInclusionLevel(0),
fOwnWs(true){
fWs = new RooWorkspace("hlfactory_ws",true);
fSigBkgPdfNames.SetOwner();
fBkgPdfNames.SetOwner();
fDatasetsNames.SetOwner();
}
HLFactory::~HLFactory(){
if (fComboSigBkgPdf!=NULL)
delete fComboSigBkgPdf;
if (fComboBkgPdf!=NULL)
delete fComboBkgPdf;
if (fComboDataset!=NULL)
delete fComboDataset;
if (fComboCat!=NULL)
delete fComboCat;
if (fOwnWs)
delete fWs;
}
int HLFactory::AddChannel(const char* label,
const char* SigBkgPdfName,
const char* BkgPdfName,
const char* DatasetName){
if (fCombinationDone){
std::cerr << "Cannot add anymore channels. "
<< "Combination already carried out.\n";
return -1;
}
if (SigBkgPdfName!=0){
if (fWs->pdf(SigBkgPdfName)==NULL){
std::cerr << "Pdf " << SigBkgPdfName << " not found in workspace!\n";
return -1;
}
TObjString* name = new TObjString(SigBkgPdfName);
fSigBkgPdfNames.Add(name);
}
if (BkgPdfName!=0){
if (fWs->pdf(BkgPdfName)==NULL){
std::cerr << "Pdf " << BkgPdfName << " not found in workspace!\n";
return -1;
}
TObjString* name = new TObjString(BkgPdfName);
fBkgPdfNames.Add(name);
}
if (DatasetName!=0){
if (fWs->data(DatasetName)==NULL){
std::cerr << "Dataset " << DatasetName << " not found in workspace!\n";
return -1;
}
TObjString* name = new TObjString(DatasetName);
fDatasetsNames.Add(name);
}
if (label!=0){
TObjString* name = new TObjString(label);
fLabelsNames.Add(name);
}
return 0;
}
RooAbsPdf* HLFactory::GetTotSigBkgPdf(){
if (fSigBkgPdfNames.GetSize()==0)
return 0;
if (fComboSigBkgPdf!=NULL)
return fComboSigBkgPdf;
if (!fNamesListsConsistent())
return NULL;
if (fSigBkgPdfNames.GetSize()==1){
TString name(((TObjString*)fSigBkgPdfNames.At(0))->String());
fComboSigBkgPdf=fWs->pdf(name);
return fComboSigBkgPdf;
}
if (!fCombinationDone)
fCreateCategory();
RooArgList pdfs("pdfs");
TIterator* it=fSigBkgPdfNames.MakeIterator();
TObjString* ostring;
TObject* obj;
while ((obj = it->Next())){
ostring=(TObjString*) obj;
pdfs.add( *(fWs->pdf(ostring->String())) );
}
delete it;
TString name(GetName());
name+="_sigbkg";
TString title(GetName());
title+="_sigbkg";
fComboSigBkgPdf=
new RooSimultaneous(name,
title,
pdfs,
*fComboCat);
return fComboSigBkgPdf;
}
RooAbsPdf* HLFactory::GetTotBkgPdf(){
if (fBkgPdfNames.GetSize()==0)
return 0;
if (fComboBkgPdf!=NULL)
return fComboBkgPdf;
if (!fNamesListsConsistent())
return NULL;
if (fBkgPdfNames.GetSize()==1){
fComboBkgPdf=fWs->pdf(((TObjString*)fBkgPdfNames.First())->String());
return fComboBkgPdf;
}
if (!fCombinationDone)
fCreateCategory();
RooArgList pdfs("pdfs");
TIterator* it = fBkgPdfNames.MakeIterator();
TObjString* ostring;
TObject* obj;
while ((obj = it->Next())){
ostring=(TObjString*) obj;
pdfs.add( *(fWs->pdf(ostring->String())) );
}
TString name(GetName());
name+="_bkg";
TString title(GetName());
title+="_bkg";
fComboBkgPdf=
new RooSimultaneous(name,
title,
pdfs,
*fComboCat);
return fComboBkgPdf;
}
RooDataSet* HLFactory::GetTotDataSet(){
if (fDatasetsNames.GetSize()==0)
return 0;
if (fComboDataset!=NULL)
return fComboDataset;
if (!fNamesListsConsistent())
return NULL;
if (fDatasetsNames.GetSize()==1){
fComboDataset=(RooDataSet*)fWs->data(((TObjString*)fDatasetsNames.First())->String());
return fComboDataset;
}
if (!fCombinationDone)
fCreateCategory();
TIterator* it = fDatasetsNames.MakeIterator();
TObjString* ostring;
TObject* obj = it->Next();
ostring = (TObjString*) obj;
fComboDataset = (RooDataSet*) fWs->data(ostring->String()) ;
if (!fComboDataset) return NULL;
fComboDataset->Print();
TString dataname(GetName());
fComboDataset = new RooDataSet(*fComboDataset,dataname+"_TotData");
int catindex=0;
fComboCat->setIndex(catindex);
fComboDataset->addColumn(*fComboCat);
while ((obj = it->Next())){
ostring=(TObjString*) obj;
catindex++;
RooDataSet * data = (RooDataSet*)fWs->data(ostring->String());
if (!data) return NULL;
RooDataSet* dummy = new RooDataSet(*data,"");
fComboCat->setIndex(catindex);
fComboCat->Print();
dummy->addColumn(*fComboCat);
fComboDataset->append(*dummy);
delete dummy;
}
delete it;
return fComboDataset;
}
RooCategory* HLFactory::GetTotCategory(){
if (fComboCat!=NULL)
return fComboCat;
if (!fNamesListsConsistent())
return NULL;
if (!fCombinationDone)
fCreateCategory();
return fComboCat;
}
int HLFactory::ProcessCard(const char* filename){
return fReadFile(filename,0);
}
int HLFactory::fReadFile(const char*fileName, bool is_included){
if (is_included)
fInclusionLevel+=1;
else
fInclusionLevel=0;
const int maxDeepness=50;
if (fInclusionLevel>maxDeepness){
TString warning("The inclusion stack is deeper than ");
warning+=maxDeepness;
warning+=". Is this a recursive inclusion?";
Warning("fReadFile", "%s", warning.Data());
}
std::ifstream ifile(fileName);
if(ifile.fail()){
TString error("File ");
error+=fileName;
error+=" could not be opened.";
Error("fReadFile", "%s", error.Data());
return -1;
}
TString ifileContent("");
ifileContent.ReadFile(ifile);
ifile.close();
TString ifileContentStripped("");
TObjArray* lines_array = ifileContent.Tokenize("\n");
TIterator* lineIt=lines_array->MakeIterator();
bool in_comment=false;
TString line;
TObject* line_o;
while((line_o=(*lineIt)())){
line = (static_cast<TObjString*>(line_o))->GetString();
if (in_comment)
if (line.EndsWith("*/")){
in_comment=false;
if (fVerbose) Info("fReadFile","Out of multiline comment ...");
continue;
}
if ((line.BeginsWith("/*") && line.EndsWith("*/")) ||
line.BeginsWith("//")){
if (fVerbose) Info("fReadFile","In single line comment ...");
continue;
}
if (line.BeginsWith("/*")){
in_comment=true;
if (fVerbose) Info("fReadFile","In multiline comment ...");
continue;
}
ifileContentStripped+=line+"\n";
}
delete lines_array;
delete lineIt;
lines_array = ifileContentStripped.Tokenize(";");
lineIt=lines_array->MakeIterator();
in_comment=false;
const int nNeutrals=2;
TString neutrals[nNeutrals]={"\t"," "};
while((line_o=(*lineIt)())){
line = (static_cast<TObjString*>(line_o))->GetString();
line.Strip(TString::kBoth,' ');
line.ReplaceAll("\n","");
if (line.BeginsWith("echo")){
line = line(5,line.Length()-1);
if (fVerbose)
std::cout << "Echoing line " << line.Data() << std::endl;
std::cout << "[" << GetName() << "] echo: "
<< line.Data() << std::endl;
continue;
}
for (int i=0;i<nNeutrals;++i)
line.ReplaceAll(neutrals[i],"");
if (fVerbose) Info("fReadFile","Reading --> %s <--", line.Data());
if (line == ""){
if (fVerbose) Info("fReadFile", "%s", "Empty line: skipping ...");
continue;
}
if (line.BeginsWith("#include")){
line.ReplaceAll("#include","");
if (fVerbose) Info("fReadFile","Reading included file...");
fReadFile(line,true);
continue;
}
if (fVerbose) Info("fReadFile","Parsing the line...");
fParseLine(line);
}
delete lineIt;
delete lines_array;
return 0;
}
void HLFactory::fCreateCategory(){
fCombinationDone=true;
TString name(GetName());
name+="_category";
TString title(GetName());
title+="_category";
fComboCat=new RooCategory(name,title);
TIterator* it=fLabelsNames.MakeIterator();
TObjString* ostring;
TObject* obj;
while ((obj = it->Next())){
ostring=(TObjString*) obj;
fComboCat->defineType(ostring->String());
}
}
bool HLFactory::fNamesListsConsistent(){
if ((fSigBkgPdfNames.GetEntries()==fBkgPdfNames.GetEntries() || fBkgPdfNames.GetEntries()==0) &&
(fSigBkgPdfNames.GetEntries()==fDatasetsNames.GetEntries() || fDatasetsNames.GetEntries()==0) &&
(fSigBkgPdfNames.GetEntries()==fLabelsNames.GetEntries() || fLabelsNames.GetEntries()==0))
return true;
else{
std::cerr << "The number of datasets and models added as channels "
<< " is not the same!\n";
return false;
}
}
int HLFactory::fParseLine(TString& line){
if (fVerbose) Info("fParseLine", "Parsing line: %s", line.Data());
TString new_line("");
const int nequals = line.CountChar('=');
if (line.Contains("::") ||
nequals==0 ||
(line.Contains("[") &&
line.Contains("]") &&
nequals>0 &&
! line.Contains("(") &&
! line.Contains(")"))) {
fWs->factory(line);
return 0;
}
if (nequals==1 ||
(nequals > 1 && line.Contains("SIMUL"))){
const int equal_index=line.First('=');
const int par_index=line.First('(');
TString o_name(line(0,equal_index));
TString o_class(line(equal_index+1,par_index-equal_index-1));
TString o_descr(line(par_index+1,line.Length()-par_index-2));
if (fVerbose) Info("fParseLine", "o_name=%s o_class=%s o_descr=%s",
o_name.Data(), o_class.Data(), o_descr.Data());
if (o_class =="import"){
TObjArray* descr_array = o_descr.Tokenize(",");
const int n_descr_parts=descr_array->GetEntries();
if (n_descr_parts<2 || n_descr_parts>3)
Error("fParseLine","Import wrong syntax: cannot process %s", o_descr.Data());
TString obj_name (static_cast<TObjString*>(descr_array->At(n_descr_parts-1))->GetString());
TString ws_name("");
TString rootfile_name (static_cast<TObjString*>(descr_array->At(0))->GetString());
TFile* ifile=TFile::Open(rootfile_name);
if (ifile==0)
return 1;
if (n_descr_parts==3){
o_descr.ReplaceAll(",",":");
fWs->import(o_descr);
}
else if(n_descr_parts==2){
if (fVerbose)
Info("fParseLine","Importing %s from %s under the name of %s",
obj_name.Data(), rootfile_name.Data(), o_name.Data());
TObject* the_obj=ifile->Get(obj_name);
fWs->import(*the_obj,o_name);
}
delete ifile;
return 0;
}
new_line=o_class+"::"+o_name+"("+o_descr+")";
if (fVerbose){
std::cout << "DEBUG: line: " << line.Data() << std::endl;
std::cout << "DEBUG: new_line: " << new_line.Data() << std::endl;
}
fWs->factory(new_line);
return 0;
}
else {
fWs->factory(line);
}
return 0;
}