/////////////////////////////////////////////////////////////////////////// // // Generic Histogram container with support for groups and filling of groups // by passing a vector of data // // Authors: // Julian Book /* TOADD: reserved words, MC signals, hist classes, MC weighting, plotting! Histograms such as THxF,D, TProfile,2D,3D and n-dimensonal objects such as (THn, THnSparse) can be added via the various template functions: - AddHistogram - AddProfile + different error option (Mean,Rms,.. see: TProfile::BuildOptions(..)) - AddSparse In addition histograms can be filled with weights. Different kind of binnings can be passed (linear, arbitrary, logarithmic) either via the PairAnalysisHelper functionalities: - PairAnalysisHelper::MakeLinBinning(nbins,low,high) - PairAnalysisHelper::MakeLogBinning(nbins,low,high) - PairAnalysisHelper::MakeArbitraryBinning("1.,4.,10.,..") The 'name' and 'title;titleX;titleY;...' of the objects and all axis labels, units and names and object names are set in a consistent way and it is ensured that they are unique. Variables are added via the enumerator of PairAnalysisVarManager. Any computation of variables (max=10) are implemented according to TFormula using the TMath functions. The arguments and options in the function DrawSame and DrawTaskSame provide various and powerful possibilities to plot and post-process histograms. Examples: // 1 dimensional AddHistogram("Event", PairAnalysisHelper::MakeLinBinning(200,-1.,1.), PairAnalysisVarManager::kXvPrim); AddHistogram("Event", PairAnalysisHelper::MakeLinBinning(200,0.0,10.0), "VtxChi/VtxNDF"); // using formulas AddHistogram("Track", PairAnalysisHelper::MakeLogBinning(400,0,4.), PairAnalysisVarManager::kPt); // logarithmic AddHistogram("Hit.TRD",PairAnalysisHelper::MakeLinBinning(105,-1.e+2,5.), "ElossdEdx-ElossMC"); // use MC truth // 2 dimensional AddHistogram("Track", PairAnalysisHelper::MakeLinBinning(100,-0.5,3.), PairAnalysisVarManager::kY, PairAnalysisHelper::MakeLinBinning(125,0,5.), PairAnalysisVarManager::kPt); AddProfile( "Track", PairAnalysisHelper::MakeLogBinning(50,0.05,10.), PairAnalysisVarManager::kP, PairAnalysisVarManager::kTRDHits, "I"); // Y-average ... */ // // /////////////////////////////////////////////////////////////////////////// #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include "PairAnalysis.h" #include "PairAnalysisHelper.h" #include "PairAnalysisSignalMC.h" #include "PairAnalysisStyler.h" #include "PairAnalysisVarManager.h" #include "PairAnalysisHistos.h" #include "PairAnalysisMetaData.h" ClassImp(PairAnalysisHistos) PairAnalysisHistos::PairAnalysisHistos() : // TCollection(), TNamed("PairAnalysisHistos","PairAnalysis Histogram Container"), fHistoList(), fMetaData(0x0), fList(0x0), fUsedVars(new TBits(PairAnalysisVarManager::kNMaxValues)), fReservedWords(new TString("Hit;Track;Pair")), fPrecision(kFloat) { // // Default constructor // fHistoList.SetOwner(kTRUE); fHistoList.SetName("PairAnalysis_Histos"); PairAnalysisVarManager::InitFormulas(); PairAnalysisStyler::LoadStyle(); } //_____________________________________________________________________________ PairAnalysisHistos::PairAnalysisHistos(const char* name, const char* title) : // TCollection(), TNamed(name, title), fHistoList(), fMetaData(0x0), fList(0x0), fUsedVars(new TBits(PairAnalysisVarManager::kNMaxValues)), fReservedWords(new TString("Hit;Track;Pair")), fPrecision(kFloat) { // // TNamed constructor // fHistoList.SetOwner(kTRUE); fHistoList.SetName(name); PairAnalysisVarManager::InitFormulas(); PairAnalysisStyler::LoadStyle(); } //_____________________________________________________________________________ PairAnalysisHistos::~PairAnalysisHistos() { // // Destructor // fHistoList.Clear(); if (fUsedVars) delete fUsedVars; if (fList) fList->Clear(); if (fMetaData) delete fMetaData; delete fReservedWords; } //_____________________________________________________________________________ void PairAnalysisHistos::UserHistogram(const char* histClass, Int_t ndim, TObjArray *limits, UInt_t *vars, UInt_t valTypeW) { // // Histogram creation n>3 dimension only with non-linear binning // Bool_t isOk=kTRUE; isOk&=(ndim<21 && ndim>3); if(!isOk) { Warning("UserHistogram","Array sizes should be between 3 and 20. Not adding Histogram to '%s'.", histClass); return; } isOk&=(ndim==limits->GetEntriesFast()); if(!isOk) return; // set automatic histo name TString name; for(Int_t iv=0; iv < ndim; iv++) name+=Form("%s_",PairAnalysisVarManager::GetValueName(vars[iv])); name.Resize(name.Length()-1); isOk&=IsHistogramOk(histClass,name); THnD *hist; Int_t bins[ndim]; if (isOk) { // get number of bins for(Int_t idim=0 ;idimAt(idim); bins[idim]=vec->GetNrows()-1; } hist=new THnD(name.Data(),"", ndim, bins, 0x0, 0x0); // set binning for(Int_t idim=0 ;idimAt(idim); hist->SetBinEdges(idim,vec->GetMatrixArray()); } // store variales in axes StoreVariables(hist, vars); hist->SetUniqueID(valTypeW); // store weighting variable // store which variables are used for(Int_t i=0; iSetBitNumber(vars[i],kTRUE); fUsedVars->SetBitNumber(valTypeW,kTRUE); Bool_t isReserved=fReservedWords->Contains(histClass); if (isReserved) UserHistogramReservedWords(histClass, hist); else UserHistogram(histClass, hist); } } //_____________________________________________________________________________ void PairAnalysisHistos::AddSparse(const char* histClass, Int_t ndim, TObjArray *limits, UInt_t *vars, UInt_t valTypeW) { // // THnSparse creation with non-linear binning // Bool_t isOk=kTRUE; isOk&=(ndim==limits->GetEntriesFast()); if(!isOk) return; // set automatic histo name TString name; for(Int_t iv=0; iv < ndim; iv++) name+=Form("%s_",PairAnalysisVarManager::GetValueName(vars[iv])); name.Resize(name.Length()-1); isOk&=IsHistogramOk(histClass,name); // THnSparseD *hist; PairAnalysisHn *hist; Int_t bins[ndim]; if (isOk) { // get number of bins for(Int_t idim=0 ;idimAt(idim); bins[idim]=vec->GetNrows()-1; } // hist=new THnSparseD(name.Data(),"", ndim, bins, 0x0, 0x0); hist=new PairAnalysisHn(name.Data(),"", ndim, bins, 0x0, 0x0); // set binning for(Int_t idim=0 ;idimAt(idim); hist->SetBinEdges(idim,vec->GetMatrixArray()); } // store variales in axes StoreVariables(hist, vars); hist->SetUniqueID(valTypeW); // store weighting variable // store which variables are used for(Int_t i=0; iSetBitNumber(vars[i],kTRUE); fUsedVars->SetBitNumber(valTypeW,kTRUE); Bool_t isReserved=fReservedWords->Contains(histClass); if (isReserved) UserHistogramReservedWords(histClass, hist); else UserHistogram(histClass, hist); } } //_____________________________________________________________________________ void PairAnalysisHistos::AddSparse(const char* histClass, Int_t ndim, TObjArray *limits, TFormula **vars, UInt_t valTypeW) { // // THnSparse creation with non-linear binning // Bool_t isOk=kTRUE; isOk&=(ndim==limits->GetEntriesFast()); if(!isOk) return; // set automatic histo name TString name; for(Int_t iv=0; iv < ndim; iv++) name+=Form("%s_",(PairAnalysisHelper::GetFormulaName(vars[iv])).Data()); name.Resize(name.Length()-1); name.ReplaceAll("f(",""); name.ReplaceAll(")",""); name.Prepend("f("); name.Append(")"); isOk&=IsHistogramOk(histClass,name); // THnSparseD *hist; PairAnalysisHn *hist; Int_t bins[ndim]; if (isOk) { // get number of bins for(Int_t idim=0 ;idimAt(idim); bins[idim]=vec->GetNrows()-1; } // hist=new THnSparseD(name.Data(),"", ndim, bins, 0x0, 0x0); hist=new PairAnalysisHn(name.Data(),"", ndim, bins, 0x0, 0x0); // set binning for(Int_t idim=0 ;idimAt(idim); hist->SetBinEdges(idim,vec->GetMatrixArray()); } // add formulas to list of functions for(Int_t idim=0 ;idimSetName(Form("axis%dFormula",idim)); hist->GetListOfFunctions()->Add( vars[idim] ); } // store variales in axes // StoreVariables(hist, vars); hist->SetUniqueID(valTypeW); // store weighting variable // adapt the name and title of the histogram in case they are empty AdaptNameTitle(hist, histClass); // store which variables are used //for(Int_t i=0; iSetBitNumber(vars[i],kTRUE); fUsedVars->SetBitNumber(valTypeW,kTRUE); Bool_t isReserved=fReservedWords->Contains(histClass); if (isReserved) UserHistogramReservedWords(histClass, hist); else UserHistogram(histClass, hist); } } //_____________________________________________________________________________ TString PairAnalysisHistos::UserHistogram(const char* histClass, TObject* hist) { // // Add any type of user histogram // //special case for the calss Pair. where histograms will be created for all pair classes Bool_t isReserved=fReservedWords->Contains(histClass); if (isReserved) { UserHistogramReservedWords(histClass, hist); return hist->GetName(); } // if (!IsHistogramOk(histClass,hist->GetName())) return hist->GetName(); // THashList *classTable=(THashList*)fHistoList.FindObject(histClass); // // hist->SetDirectory(0); // store variables axis UInt_t valType[20] = {0}; // extract variables from axis FillVarArray(hist, valType); // TODO: implement THn/PairAnalysisHn // change to mc truth variables if available TString hclass=histClass; if(hclass.Contains("MCtruth") && !(hist->IsA()==PairAnalysisHn::Class()) ) { for(Int_t i=0;i<2;i++) { /// protection for changes of variable enum (in the postprocessing) --> use axistname indentification // UInt_t valTypeFromTitle = 0; if(!i) valType[i] = PairAnalysisVarManager::GetValueType(((TH1*)hist)->GetXaxis()->GetName()); else valType[i] = PairAnalysisVarManager::GetValueType(((TH1*)hist)->GetYaxis()->GetName()); //Printf("SWITCH TO MC: before: %d %s, (via axis name %d) ---->",valType[i],PairAnalysisVarManager::GetValueName(valType[i]),valTypeFromTitle); valType[i] = PairAnalysisVarManager::GetValueTypeMC(valType[i]); // if theres no corresponding MCtruth variable, skip adding this histogram // if(valType[i] < PairAnalysisVarManager::kNMaxValues && valType[i]>0) return hist->GetName(); if(valType[i] < PairAnalysisVarManager::kPairMax && valType[i]>0) return hist->GetName(); // request filling of mc variable fUsedVars->SetBitNumber(valType[i],kTRUE); // Printf("after: %d %s",valType[i],PairAnalysisVarManager::GetValueName(valType[i])); } StoreVariables(hist, valType); hist->SetUniqueID(valType[19]); // store weighting variable // check for formulas if(hist->InheritsFrom(TH1::Class())) { TIter next(((TH1*)hist)->GetListOfFunctions()); TFormula *f=0; while ( (f=dynamic_cast(next()) ) ) { for(Int_t i=0; iGetNpar(); i++) { Int_t parMC=PairAnalysisVarManager::GetValueTypeMC(f->GetParameter(i)); // if theres none corresponding MCtruth variable, skip adding this histogram if(parMC < PairAnalysisVarManager::kNMaxValues) return hist->GetName(); f->SetParameter( i, parMC ); f->SetParName( i, PairAnalysisVarManager::GetValueName(parMC) ); fUsedVars->SetBitNumber(parMC,kTRUE); } } // change histogram key according to mctruth information AdaptNameTitle((TH1*)hist, histClass); } } else { StoreVariables(hist, valType); hist->SetUniqueID(valType[19]); // store weighting variable } // add histogram to class if (!IsHistogramOk(histClass,hist->GetName())) return hist->GetName(); THashList *classTable=(THashList*)fHistoList.FindObject(histClass); // hist->SetDirectory(0); if(classTable) classTable->Add(hist); return hist->GetName(); } //_____________________________________________________________________________ TH1 *PairAnalysisHistos::GetTHist(const char* histClass, const char *name, const char* title, const TVectorD * const binsX, const TVectorD * const binsY, const TVectorD * const binsZ ) { // // retrieve n-dimensional Hist depending on arguments // Bool_t isOk=kTRUE; isOk&=IsHistogramOk(histClass,name); isOk&=(binsX!=0x0); if (!isOk) return 0x0; switch(fPrecision) { case kFloat: if(!binsY) return (new TH1F(name,title,binsX->GetNrows()-1,binsX->GetMatrixArray()) ); else if(!binsZ) return (new TH2F(name,title,binsX->GetNrows()-1,binsX->GetMatrixArray(), binsY->GetNrows()-1,binsY->GetMatrixArray()) ); else return (new TH3F(name,title,binsX->GetNrows()-1,binsX->GetMatrixArray(), binsY->GetNrows()-1,binsY->GetMatrixArray(), binsZ->GetNrows()-1,binsZ->GetMatrixArray()) ); break; case kDouble: if(!binsY) return (new TH1D(name,title,binsX->GetNrows()-1,binsX->GetMatrixArray()) ); else if(!binsZ) return (new TH2D(name,title,binsX->GetNrows()-1,binsX->GetMatrixArray(), binsY->GetNrows()-1,binsY->GetMatrixArray()) ); else return (new TH3D(name,title,binsX->GetNrows()-1,binsX->GetMatrixArray(), binsY->GetNrows()-1,binsY->GetMatrixArray(), binsZ->GetNrows()-1,binsZ->GetMatrixArray()) ); break; default: return 0x0; break; } } //_____________________________________________________________________________ TH1 *PairAnalysisHistos::GetTProf(const char* histClass, const char *name, const char* title, const TVectorD * const binsX, const TVectorD * const binsY, const TVectorD * const binsZ, TString option) { // // retrieve n-dimensional profile histogram with error options depending on arguments // Bool_t isOk=kTRUE; isOk&=IsHistogramOk(histClass,name); isOk&=(binsX!=0x0); if (!isOk) return 0x0; // parse error option TString opt=""; Double_t pmin=0., pmax=0.; if(!option.IsNull()) { TObjArray *arr=option.Tokenize(";"); arr->SetOwner(); opt=((TObjString*)arr->At(0))->GetString(); if(arr->GetEntriesFast()>1) pmin=(((TObjString*)arr->At(1))->GetString()).Atof(); if(arr->GetEntriesFast()>2) pmax=(((TObjString*)arr->At(2))->GetString()).Atof(); delete arr; } // build profile with error options and return it TH1 *prof=0x0; if(!binsY) return (new TProfile(name,title,binsX->GetNrows()-1,binsX->GetMatrixArray(),pmin,pmax,opt.Data()) ); else if(!binsZ) { prof = new TProfile2D(name,title,binsX->GetNrows()-1,binsX->GetMatrixArray(), binsY->GetNrows()-1,binsY->GetMatrixArray()); ((TProfile2D*)prof)->BuildOptions(pmin,pmax,opt.Data()); return prof; } else { prof = new TProfile3D(name,title,binsX->GetNrows()-1,binsX->GetMatrixArray(), binsY->GetNrows()-1,binsY->GetMatrixArray(), binsZ->GetNrows()-1,binsZ->GetMatrixArray()); ((TProfile3D*)prof)->BuildOptions(pmin,pmax,opt.Data()); return prof; } } //_____________________________________________________________________________ TFormula *PairAnalysisHistos::GetFormula(const char *name, const char* formula) { // // build a TFormula object // TFormula *form = new TFormula(name,formula); // compile function if(form->Compile()) return 0x0; //set parameter/variable identifier for(Int_t i=0; iGetNpar(); i++) { form->SetParName( i, PairAnalysisVarManager::GetValueName(form->GetParameter(i)) ); fUsedVars->SetBitNumber((Int_t)form->GetParameter(i),kTRUE); } return form; } //_____________________________________________________________________________ void PairAnalysisHistos::AddClass(const char* histClass) { // // Add a class of histograms // Several classes can be added by separating them by a ';' e.g. 'class1;class2;class3' // TString hists(histClass); TObjArray *arr=hists.Tokenize(";"); TIter next(arr); TObject *o=0; while ( (o=next()) ){ if (fHistoList.FindObject(o->GetName())){ Warning("AddClass","Cannot create class '%s' it already exists.",histClass); continue; } if (fReservedWords->Contains(o->GetName())){ Error("AddClass","Pair is a reserved word, please use another name"); continue; } THashList *table=new THashList; table->SetOwner(kTRUE); table->SetName(o->GetName()); fHistoList.Add(table); } delete arr; } //_____________________________________________________________________________ void PairAnalysisHistos::FillClass(TString histClass, const Double_t *values) { // // Fill class 'histClass' (by name) // THashList *classTable=(THashList*)fHistoList.FindObject(histClass.Data()); if (!classTable){ // Warning("FillClass","Cannot fill class '%s' its not defined.",histClass.Data()); return; } TIter nextHist(classTable); TObject *obj=0; while ( (obj=(TObject*)nextHist()) ) FillValues(obj, values); return; } //_____________________________________________________________________________ void PairAnalysisHistos::UserHistogramReservedWords(const char* histClass, const TObject *hist) { // // Creation of histogram for all pair or track types // TString title(hist->GetTitle()); TIter nextClass(&fHistoList); THashList *l=0; while ( (l=static_cast(nextClass())) ){ TString name(l->GetName()); if (name.Contains(histClass)){ TObject *h=hist->Clone(); // Tobject has no function SetDirectory, didn't we need this??? // h->SetDirectory(0); if(h->InheritsFrom(TH1::Class())) ((TH1*)h)->SetTitle(Form("%s %s",title.Data(),l->GetName())); else ((THnBase*)h)->SetTitle(Form("%s %s",title.Data(),l->GetName())); UserHistogram(l->GetName(),h); } } delete hist; return; } //_____________________________________________________________________________ void PairAnalysisHistos::DumpToFile(const char* file) { // // Dump the histogram list to a newly created root file // TFile f(file,"recreate"); fHistoList.Write(fHistoList.GetName(),TObject::kSingleKey); f.Close(); } //_____________________________________________________________________________ TObject* PairAnalysisHistos::GetHist(const char* histClass, const char* name) const { // // return object 'name' in 'histClass' // THashList *classTable=(THashList*)fHistoList.FindObject(histClass); if (!classTable) return 0x0; return classTable->FindObject(name); } //_____________________________________________________________________________ TH1* PairAnalysisHistos::GetHistogram(const char* histClass, const char* name) const { // // return histogram 'name' in 'histClass' // return ((TH1*) GetHist(histClass, name)); } //_____________________________________________________________________________ TObject* PairAnalysisHistos::GetHist(const char* cutClass, const char* histClass, const char* name) const { // // return object from list of list of histograms // this function is thought for retrieving histograms if a list of PairAnalysisHistos is set // if (!fList) return 0x0; THashList *h=dynamic_cast(fList->FindObject(cutClass)); if (!h)return 0x0; THashList *classTable=dynamic_cast(h->FindObject(histClass)); if (!classTable) return 0x0; return classTable->FindObject(name); } //_____________________________________________________________________________ TH1* PairAnalysisHistos::GetHistogram(const char* cutClass, const char* histClass, const char* name) const { // // return histogram from list of list of histograms // this function is thought for retrieving histograms if a list of PairAnalysisHistos is set // return ((TH1*) GetHist(cutClass, histClass, name)); } //_____________________________________________________________________________ void PairAnalysisHistos::Draw(const Option_t* option) { // // Draw histograms // TString drawStr(option); TObjArray *arr=drawStr.Tokenize(";"); arr->SetOwner(); TIter nextOpt(arr); TString drawClasses; TObjString *ostr=0x0; TString currentOpt; TString testOpt; while ( (ostr=(TObjString*)nextOpt()) ){ currentOpt=ostr->GetString(); currentOpt.Remove(TString::kBoth,'\t'); currentOpt.Remove(TString::kBoth,' '); testOpt="classes="; if ( currentOpt.Contains(testOpt.Data()) ){ drawClasses=currentOpt(testOpt.Length(),currentOpt.Length()); } } delete arr; drawStr.ToLower(); //optionsfList // Bool_t same=drawOpt.Contains("same"); //FIXME not yet implemented TCanvas *c=0x0; if (gVirtualPS) { if (!gPad){ Error("Draw","When writing to a file you have to create a canvas before opening the file!!!"); return; } c=gPad->GetCanvas(); c->cd(); // c=new TCanvas; } TIter nextClass(&fHistoList); THashList *classTable=0; // Bool_t first=kTRUE; while ( (classTable=(THashList*)nextClass()) ){ //test classes option if (!drawClasses.IsNull() && !drawClasses.Contains(classTable->GetName())) continue; //optimised division Int_t nPads = classTable->GetEntries(); Int_t nCols = (Int_t)TMath::Ceil( TMath::Sqrt(nPads) ); Int_t nRows = (Int_t)TMath::Ceil( (Double_t)nPads/(Double_t)nCols ); //create canvas if (!gVirtualPS){ TString canvasName; canvasName.Form("c%s_%s",GetName(),classTable->GetName()); c=(TCanvas*)gROOT->FindObject(canvasName.Data()); if (!c) c=new TCanvas(canvasName.Data(),Form("%s: %s",GetName(),classTable->GetName())); c->Clear(); } else { // if (first){ // first=kFALSE; // if (nPads>1) gVirtualPS->NewPage(); // } else { if (nPads>1) c->Clear(); // } } if (nCols>1||nRows>1) c->Divide(nCols,nRows); //loop over histograms and draw them TIter nextHist(classTable); Int_t iPad=0; TH1 *h=0; while ( (h=(TH1*)nextHist()) ){ TString drawOpt; if ( (h->InheritsFrom(TH2::Class())) ) drawOpt="colz"; if (nCols>1||nRows>1) c->cd(++iPad); if ( TMath::Abs(h->GetXaxis()->GetBinWidth(1)-h->GetXaxis()->GetBinWidth(2))>1e-10 ) gPad->SetLogx(); if ( TMath::Abs(h->GetYaxis()->GetBinWidth(1)-h->GetYaxis()->GetBinWidth(2))>1e-10 ) gPad->SetLogy(); if ( TMath::Abs(h->GetZaxis()->GetBinWidth(1)-h->GetZaxis()->GetBinWidth(2))>1e-10 ) gPad->SetLogz(); TString histOpt=h->GetOption(); histOpt.ToLower(); if (histOpt.Contains("logx")) gPad->SetLogx(); if (histOpt.Contains("logy")) gPad->SetLogy(); if (histOpt.Contains("logz")) gPad->SetLogz(); histOpt.ReplaceAll("logx",""); histOpt.ReplaceAll("logy",""); histOpt.ReplaceAll("logz",""); h->Draw(drawOpt.Data()); } if (gVirtualPS) { c->Update(); } } // if (gVirtualPS) delete c; } //_____________________________________________________________________________ void PairAnalysisHistos::Print(const Option_t* option) const { // // Print classes and histograms // TString optString(option); if (optString.IsNull()) PrintStructure(); } //_____________________________________________________________________________ void PairAnalysisHistos::PrintStructure() const { // // Print classes and histograms in the class to stdout // if (!fList){ TIter nextClass(&fHistoList); THashList *classTable=0; while ( (classTable=(THashList*)nextClass()) ){ TIter nextHist(classTable); TObject *o=0; Printf("+ %s\n",classTable->GetName()); while ( (o=nextHist()) ) Printf("| ->%s\n",o->GetName()); } } else { TIter nextCutClass(fList); THashList *cutClass=0x0; while ( (cutClass=(THashList*)nextCutClass()) ) { TString cla=cutClass->GetName(); if(cla.Contains("QAcuts")) continue; Printf("+ %s\n",cutClass->GetName()); TIter nextClass(cutClass); THashList *classTable=0; while ( (classTable=(THashList*)nextClass()) ){ TIter nextHist(classTable); TObject *o=0; Printf("| + %s\n",classTable->GetName()); while ( (o=nextHist()) ) Printf("| | ->%s\n",o->GetName()); } } } } //_____________________________________________________________________________ void PairAnalysisHistos::SetHistogramList(THashList &list, Bool_t setOwner/*=kTRUE*/) { // // set histogram classes and histograms to this instance. It will take onwnership! // ResetHistogramList(); TString name(GetName()); // if (name == "PairAnalysisHistos") SetName(list.GetName()); TIter next(&list); TObject *o; while ( (o=next()) ){ fHistoList.Add(o); } if (setOwner){ list.SetOwner(kFALSE); fHistoList.SetOwner(kTRUE); fHistoList.SetName(list.GetName()); } else { fHistoList.SetOwner(kFALSE); fHistoList.SetName(list.GetName()); } } //_____________________________________________________________________________ Bool_t PairAnalysisHistos::SetCutClass(const char* cutClass) { // // Assign histogram list according to cutClass // if (!fList) return kFALSE; ResetHistogramList(); THashList *h=dynamic_cast(fList->FindObject(cutClass)); if (!h) { Warning("SetCutClass","cutClass '%s' not found", cutClass); return kFALSE; } SetHistogramList(*h,kFALSE); return kTRUE; } //_____________________________________________________________________________ Bool_t PairAnalysisHistos::IsHistogramOk(const char* histClass, const char* name) { // // check whether the histogram class exists and the histogram itself does not exist yet // Bool_t isReserved=fReservedWords->Contains(histClass); if (!fHistoList.FindObject(histClass)&&!isReserved){ Warning("IsHistogramOk","Cannot create histogram. Class '%s' not defined. Please create it using AddClass before.",histClass); return kFALSE; } if (GetHist(histClass,name)){ Warning("IsHistogramOk","Cannot create histogram '%s' in class '%s': It already exists!",name,histClass); return kFALSE; } return kTRUE; } // //_____________________________________________________________________________ // TIterator* PairAnalysisHistos::MakeIterator(Bool_t dir) const // { // // // // // // // return new TListIter(&fHistoList, dir); // } //_____________________________________________________________________________ void PairAnalysisHistos::ReadFromFile(const char* file, const char *task, const char *config) { // // Read histos from file // TFile f(file); TIter nextKey(f.GetListOfKeys()); TKey *key=0; while ( (key=(TKey*)nextKey()) ){ TString name=key->GetName(); // check for meta data if(name.Contains(Form("PairAnalysisMetaData_%s",task))) { fMetaData=new PairAnalysisMetaData(); fMetaData->SetMetaData(*dynamic_cast(f.Get(key->GetName())),kFALSE); } // check for histos if(!name.Contains(Form("PairAnalysisHistos_%s",task))) continue; if(!strlen(task) && !name.Contains(task)) continue; TObject *o=f.Get(key->GetName()); TList *list=dynamic_cast(o); if (!list) continue; else fList=list; THashList *listCfg=dynamic_cast(list->FindObject(config)); if (!listCfg) continue; SetHistogramList(*listCfg); break; } f.Close(); // load style PairAnalysisStyler::LoadStyle(); } //_____________________________________________________________________________ TObjArray* PairAnalysisHistos::DrawTaskSame(TString histName, TString opt, TString histClassDenom, TString taskDenom) { /// /// Draw histograms of different tasks into the same canvas /// /// additional plotting options to DrawSame: /// "selcfg": option to (de-)select certain task/config specified in "histClassDenom", /// delimiters ':;,' are recognized, exclusion done by prepending '!' /// "cutstep": histograms of the cut steps are compared for a single task (see AnalysisFilter::AddHistos) /// /// if argument "taskDenom" is provided histogram ratios to this task/config are calculated TObjArray *selections = histClassDenom.Tokenize(":;,"); opt.ToLower(); TString optString(opt); Bool_t optGoff =optString.Contains("goff"); Bool_t optEff = optString.Contains("eff"); Bool_t optCutStep = optString.Contains("cutstep"); Bool_t optSelCfg = optString.Contains("selcfg"); opt.ReplaceAll("selcfg",""); fHistoList.SetOwner(kFALSE); /// output array TObjArray *arr=NULL; if(optGoff) Info("DrawTaskSame","graphics option off, collect an array"); { // arr=(TObjArray*)gROOT->FindObject(Form("%s",histName.Data())); //arr=(TObjArray*)gROOT->FindObject(GetName()); // if(arr) { Warning("DrawTaskSame","Clear existing array %s",GetName()); arr->ls(); arr->Clear(); } // else arr = new TObjArray(); // arr->SetName(Form("%s",histName.Data())); arr->SetName(GetName()); arr->SetOwner(kFALSE); } /// legend TString legendname = GetName(); // printf("DrawTaskSame: legendname %s\n",legendname.Data()); /// reserved words TObjArray *reservedWords = fReservedWords->Tokenize(":;"); /// find cut step lists if option 'cutstep' is active if(optCutStep) { TString cutstepTask = fHistoList.GetName(); THashList *listCutStep=dynamic_cast(fList->FindObject(cutstepTask)); if(listCutStep) fList=listCutStep; } /// find denominator list 'taskDenom' THashList *listDenom = dynamic_cast(fList->FindObject(taskDenom.Data())); if(listDenom) opt+="div"; /// iterate over all tasks and their lists in the file TIter nextCfg(fList); THashList *listCfg=0; while ( (listCfg=static_cast(nextCfg())) ){ TString lname=listCfg->GetName(); /// exclude QAcuts if( lname.Contains("QAcuts_") ) continue; Info("DrawTaskSame"," Task name %s ",lname.Data()); /// skip same config if argument 'taskDenom' is set if(!optEff && listDenom && lname.EqualTo(taskDenom)) continue; /// skip configs according to selection done via histClassDenom and option 'selcfg' if(optSelCfg && selections->GetEntriesFast()) { Bool_t pass = kFALSE; for(Int_t is=0; isGetEntriesFast(); is++) { Bool_t testIgnore = kFALSE; TString raw = ((TObjString*)selections->At(is))->GetString(); // printf("sel: '%s' \n", histClassDenom.Data()); /// check if selection string contains reserved word /// if true ignore it for config selection for(Int_t ir=0; irGetEntriesFast(); ir++) { testIgnore = raw.Contains( ((TObjString*)reservedWords->At(ir))->GetString()); if(testIgnore) break; } if(testIgnore) continue; TString srch = raw; srch.ReplaceAll("!",""); Bool_t optExclSel = !(srch.EqualTo(raw)); // exclude or not // printf("selection %d/%d: raw '%s' \t search: '%s' exclude: %d compare to: '%s' \t pass %d \n", // is,selections->GetEntriesFast()-1,raw.Data(),srch.Data(),optExclSel,lname.Data(),pass); // printf("\t decision !(%d^%d)=%d \n",(!lname.EqualTo(srch)),(optExclSel),!(!lname.EqualTo(srch))^(optExclSel)); /// decision if( !(!lname.EqualTo(srch))^(optExclSel) ) pass = kTRUE; /// same config found - direct decision (skip other selections) if( lname.EqualTo(srch) ) { pass = !(!lname.EqualTo(srch))^(optExclSel); // printf("\t same config found -> stop with pass: %d \n",pass); break; } // remove config selection string abgearbeitet // histClassDenom.ReplaceAll(raw.Data(),""); } if(!pass) continue; } /// keep only legend name!=task name if(lname.EqualTo(legendname)) legendname=""; /// update current histogram list ResetHistogramList(); TIter next(listCfg); Int_t idx=0; TObject *o; while ( (o=next()) ){ fHistoList.AddAt(o,idx++); } //fHistoList.Print(); /// adapt name for legend SetName(listCfg->GetName()); arr->AddAll( DrawSame(histName,(opt+"task").Data(),histClassDenom, listDenom) ); } /// set legend name if(optString.Contains("leg")) { // printf("DrawTaskSame: SET legendname %s\n",legendname.Data()); legendname.ReplaceAll(".",""); TList *prim = gPad->GetListOfPrimitives(); TLegend *leg=(TLegend*)prim->FindObject("TPave"); // remove legend header from legend entry TList *llist = leg->GetListOfPrimitives(); Int_t nent = llist->GetEntries(); for(Int_t il=0; il(llist->At(il)); TString lst(lent->GetLabel()); lst.ReplaceAll(legendname.Data(),""); lent->SetLabel(lst.Data()); } if(!legendname.EqualTo("none")) leg->SetHeader(""); else leg->SetHeader(legendname.Data()); PairAnalysisStyler::SetLegendAttributes(leg); // coordinates, margins, fillstyle, fontsize leg->Draw(); } /// cleanup if(selections) delete selections; if(reservedWords) delete reservedWords; if(!optGoff) delete arr; return arr; } //_____________________________________________________________________________ TObjArray* PairAnalysisHistos::DrawSame(TString histName, TString option, TString histClassDenom, THashList *listDenom) { /// /// Draw histograms with the same histName into a canvas /// /// additional plotting options to TH1::Draw(): /// "back": iterate backwards over the histClasses /// "goff": graphics are switched off & an array with the objects is returned /// "task": histograms of different tasks are compared (see DrawTaskSame) /// /// "div": histograms of different tasks are divided (see DrawTaskSame) INTERNAL? /// "eff": efficiencies are calculated /// "ratio": ratios of any histclass to "histClassDenom" are calculated /// "oneOver": the histogram is rescaled by 1./(bin contents) /// "stack": build and plot a stack histogram /// /// "noMc": no mc signals are plotted /// "noMcTrue": no mc truth signals are plotted /// "onlyMc": only mc signals are plotted /// "onlyMcTrue": only mc truth signals are plotted /// "sel": option to (de-)select certain histClasses specified in "histClassDenom", /// delimiters ":;," are recognized, exclusion done by prepending "!" /// /// "can": a new canvas is created with name: "c"+"histName" /// "logx,y,z": the axis are plotted in log-scale (labels are added automatically according to the range) /// "meta" adds default meta data to the pad (see PairAnalysisMetaData::DrawSame) /// "legF": (option "F":filled) legend will be created with caption=className , /// can be modified by PairAnalysisHistos::SetName("mycaption"), /// change of legend position: see PairAnalysisStyler::SetLegendAlign /// "meanX,Y": quote the mean in the legend (Y only works for TProfile a.t.m.) /// "rmsX,Y": quote the rms in the legend (Y only works for TProfile a.t.m.) /// "det": adds the detector name to the legend for detector specific histograms (e.g. hit histograms) /// /// "width": scale histograms by 1./bin width according to TH1::Scale /// "rebinstat0.XY": the objects are rebinned to have <"0.XY" stat uncertainty in each bin (max two digits) /// "rebinX": rebins the histogram along x-axis by factor "X" (for X<10) /// "sclmax": scale the histogram by 1./maximum /// "norm": histogram is normalized to 1. /// "normY": 2D-histograms are normalized in x-axis slices to 1. /// "cum": calculate cumulated sum of bins (+under- & overflow) for 1D-histograms (left-to-right) /// "cumR": calculate cumulated sum of bins (+under- & overflow)for 1D-histograms ("R"=reverse: right-to-left) /// "events": use number of selected (after cuts) events in meta data to normalize the histograms /// "smoothX": smooths the histogram along x-axis by factor "X" (for X<10) /// /// "rststy": reset style index to zero (usefull when the same color code should be plotted on top of each other) /// "geant": translate geantId to geant process names (see PairAnalysisHelper::SetGEANTBinLabels) /// "pdg": translate bin low edges into pdg label /// "pdgc": translate bin low edges into pdg label for bins with entries TString optString(option); optString.ToLower(); printf("Plot hist: '%s' class-denom/sel: '%s' \t listDenom: '%s' \t options: '%s' \n", histName.Data(), histClassDenom.Data(), (listDenom?listDenom->GetName():""), optString.Data()); Bool_t optBack =optString.Contains("back"); optString.ReplaceAll("back",""); Bool_t optGoff =optString.Contains("goff"); optString.ReplaceAll("goff",""); Bool_t optTask =optString.Contains("task"); optString.ReplaceAll("task",""); Bool_t optCutStep = optString.Contains("cutstep"); optString.ReplaceAll("cutstep",""); /// options - calulation Bool_t optDiv =optString.Contains("div"); optString.ReplaceAll("div",""); Bool_t optEff =optString.Contains("eff"); optString.ReplaceAll("eff",""); Bool_t optRatio =optString.Contains("ratio"); optString.ReplaceAll("ratio",""); Bool_t optOneOver =optString.Contains("oneover"); optString.ReplaceAll("oneover",""); /// options - selection Bool_t optNoMCtrue =optString.Contains("nomctrue"); optString.ReplaceAll("nomctrue",""); Bool_t optNoMC =optString.Contains("nomc"); optString.ReplaceAll("nomc",""); Bool_t optOnlyMCtrue=optString.Contains("onlymctrue"); optString.ReplaceAll("onlymctrue",""); Bool_t optOnlyMC =optString.Contains("onlymc"); optString.ReplaceAll("onlymc",""); Bool_t optSel =optString.Contains("sel"); optString.ReplaceAll("sel",""); /// options - representation Bool_t optCan =optString.Contains("can"); optString.ReplaceAll("can",""); Bool_t optLegFull =optString.Contains("legf"); optString.ReplaceAll("legf",""); Bool_t optLeg =optString.Contains("leg"); optString.ReplaceAll("leg",""); Bool_t optDet =optString.Contains("det"); optString.ReplaceAll("det",""); Bool_t optMeta =optString.Contains("meta"); optString.ReplaceAll("meta",""); Bool_t optWdth =optString.Contains("width"); optString.ReplaceAll("width",""); Bool_t optRbnStat =optString.Contains("rebinstat"); optString.ReplaceAll("rebinstat","stat"); Bool_t optRbn =optString.Contains("rebin"); Bool_t optSmooth =optString.Contains("smooth"); Bool_t optSclMax =optString.Contains("sclmax"); optString.ReplaceAll("sclmax",""); Bool_t optNormY =optString.Contains("normy"); optString.ReplaceAll("normy",""); Bool_t optNorm =optString.Contains("norm"); optString.ReplaceAll("norm",""); Bool_t optCumR =optString.Contains("cumr"); optString.ReplaceAll("cumr",""); Bool_t optCum =optString.Contains("cum"); optString.ReplaceAll("cum",""); Bool_t optEvt =optString.Contains("events"); optString.ReplaceAll("events",""); Bool_t optStack =optString.Contains("stack"); optString.ReplaceAll("stack",""); Bool_t optRstSty =optString.Contains("rststy"); optString.ReplaceAll("rststy",""); Bool_t optSlicesY =optString.Contains("slicesy"); optString.ReplaceAll("slicesy",""); /// options - information Bool_t optMeanX =optString.Contains("meanx"); optString.ReplaceAll("meanx",""); Bool_t optRmsX =optString.Contains("rmsx"); optString.ReplaceAll("rmsx",""); Bool_t optMeanY =optString.Contains("meany"); optString.ReplaceAll("meany",""); Bool_t optRmsY =optString.Contains("rmsy"); optString.ReplaceAll("rmsy",""); Bool_t optGeant =optString.Contains("geant"); optString.ReplaceAll("geant",""); Bool_t optPdgClean =optString.Contains("pdgc"); optString.ReplaceAll("pdgc",""); Bool_t optPdg =optString.Contains("pdg"); optString.ReplaceAll("pdg",""); /// set rebinning Int_t rbn = 2; if(optRbn) { TString rebin(optString(optString.Index("rebin",5,0,TString::kExact)+5)); if(rebin.IsDigit()) { rbn = rebin.Atoi(); optString.ReplaceAll(rebin,""); } optString.ReplaceAll("rebin",""); } /// set stat. uncertainty limit for rebinning /// NOTE: first check for one digit, then 2 digits TArrayD* limits=NULL; Double_t stat = 0.8; if(optRbnStat) { TString rebin(optString(optString.Index("stat",4,0,TString::kExact)+4,4)); if(!rebin.IsFloat()) rebin=optString(optString.Index("stat",4,0,TString::kExact)+4,3); if(rebin.IsFloat()) { stat = rebin.Atof(); // printf("rebinstat string: '%s' -> %f\n",rebin.Data(),stat); optString.ReplaceAll(rebin,""); } optString.ReplaceAll("stat",""); } /// set smoothing Int_t smth = 1; if(optSmooth) { TString smooth(optString(optString.Index("smooth",6,0,TString::kExact)+6)); if(smooth.IsDigit()) { smth = smooth.Atoi(); optString.ReplaceAll(smooth,""); } optString.ReplaceAll("smooth",""); } /// activate std option for legend and conncetions if(optLegFull) optLeg=kTRUE; if(optCumR) optCum=kTRUE; /// selection string TString select(""); TObjArray *selections = (optSel ? histClassDenom.Tokenize(":;,") : NULL); if(optSel) histClassDenom=""; /// output array TObjArray *arr=0x0; if(optGoff) { Info("DrawSame","graphics option off, collect an array"); // arr=(TObjArray*)gROOT->FindObject(Form("%s",histName.Data())); arr=(TObjArray*)gROOT->FindObject(GetName()); if(arr) { Warning("DrawSame","Clear existing array %s",GetName()); arr->ls(); arr->Clear(); } else arr = new TObjArray(); // arr->SetName(Form("%s",histName.Data())); arr->SetName(GetName()); arr->SetOwner(kFALSE); } /// stack option THStack *hs = NULL; // if(optStack) { // hs = new THStack("hs","stacked histograms"); // } /// add canvas TCanvas *c=0; if (optCan){ c=(TCanvas*)gROOT->FindObject(Form("c%s",histName.Data())); if (!c) { c=new TCanvas(Form("c%s",histName.Data()),Form("All '%s' histograms",histName.Data())); } // c->Clear(); c->cd(); } /// count number of drawn objects in pad TH1 *hFirst=0x0; TObject *obj; Int_t nobj=0; TList *prim = (gPad ? gPad->GetListOfPrimitives() : 0x0); // if(prim) prim->ls(); // if(prim && prim->GetSize()>1 && !optGoff) prim->RemoveLast(); // remove redraw axis histogram for(Int_t io=0; io< (prim ? prim->GetSize() : 0); io++) { obj=prim->At(io); if(obj->InheritsFrom(TH1::Class()) && obj!=prim->At(io+1)) nobj++; if(obj->InheritsFrom(THStack::Class()) && obj!=prim->At(io+1)) nobj++; if(nobj==1) hFirst=(TH1*)obj; } /// add or get legend TLegend *leg=0; if( (optLeg && !nobj) ) { leg=new TLegend(0. + gPad->GetLeftMargin() + gStyle->GetTickLength("Y"), 0. + gPad->GetBottomMargin()+ gStyle->GetTickLength("X"), 1. - gPad->GetRightMargin() + gStyle->GetTickLength("Y"), 1. - gPad->GetTopMargin() + gStyle->GetTickLength("X"), GetName(),"nbNDC"); if(optTask && !optCutStep) leg->SetHeader(""); } else if(nobj) { leg=(TLegend*)prim->FindObject("TPave"); // leg->SetX1(0. + gPad->GetLeftMargin() + gStyle->GetTickLength("Y")); // leg->SetY1(0. + gPad->GetBottomMargin()+ gStyle->GetTickLength("X")); // leg->SetX2(1. - gPad->GetRightMargin() + gStyle->GetTickLength("Y")); // leg->SetY2(1. - gPad->GetTopMargin() + gStyle->GetTickLength("X")); // leg->SetOption("nbNDC"); } Info("DrawSame","Basics: nobj: %d \t leg: %p",nobj,leg); /// logaritmic axes if(optString.Contains("logx")) gPad->SetLogx(); if(optString.Contains("logy")) gPad->SetLogy(); if(optString.Contains("logz")) gPad->SetLogz(); optString.ReplaceAll("logx",""); optString.ReplaceAll("logy",""); optString.ReplaceAll("logz",""); /// meta data Int_t events = 1; if(fMetaData && optEvt) fMetaData->GetMeta("events",&events); Int_t i=nobj; if(optTask && nobj) i=nobj; TListIter next(&fHistoList,(optBack?kIterBackward:kIterForward)); //TIter next(&fHistoList); THashList *classTable=0; // TH1 *hFirst=0x0; /// iterate over all histogram classes (events,pairs,tracks,hits) while ( (classTable=(THashList*)next()) ){ TString histClass=classTable->GetName(); Int_t ndel = histClass.CountChar('_'); Info("DrawSame","Search for hist: '%s' class-denom: '%s' select: '%s' \t ndel: %d \t for class: '%s'", histName.Data(), histClassDenom.Data(), select.Data(), ndel, histClass.Data() ); /// check selections done via MC options 'onlyMC', 'noMC', 'noMCtrue', 'eff' if( (optNoMC && ndel>0) || (optEff && ndel<1) || (optNoMCtrue && histClass.Contains("_MCtruth")) || (optOnlyMC && ndel<1) || (optOnlyMCtrue && ndel<2) ) continue; /// histclass selections done via option 'sel' and argument 'histClassDenom' -> into array selections if(optSel && selections->GetEntriesFast()) { Bool_t pass = kFALSE; for(Int_t is=0; isGetEntriesFast(); is++) { TString raw = ((TObjString*)selections->At(is))->GetString(); TString srch = raw; srch.ReplaceAll("!",""); Bool_t optExclSel = !(srch.EqualTo(raw)); // exclude or not /// decision if( !(!histClass.Contains(srch,TString::kIgnoreCase))^(optExclSel) ) pass = kTRUE; /// exact string found - direct decision (ignore other selections) if( histClass.EqualTo(srch,TString::kIgnoreCase) ) { pass = !(!histClass.EqualTo(srch,TString::kIgnoreCase))^(optExclSel); break; } } if(!pass) continue; } /// find the histogram 'histName' in the class table /// or the corressponding MC histogram TH1 *h=(TH1*)classTable->FindObject(histName.Data()); if (!h && !optNoMC && !optNoMCtrue) { /// TODO: improve for tprofiles, formulas and more dimensional histos TString histdenomMC = histName+"MC"; h = (TH1*) classTable->FindObject( histdenomMC.Data() ); } if (!h) continue; if (h->GetEntries()<1.) continue; /// get histClassDenom for efficiency caluclation, e.g. the MCtruth (denominator) if(optEff && !histClass.Contains("_MCtruth")) histClassDenom = histClass + "_MCtruth"; Info("DrawSame"," Hist found in histClass '%s' (search for denom '%s') ",histClass.Data(),histClassDenom.Data()); /// check if 'ratio' or 'eff' should be build if( (optEff || optRatio) /*&& !optTask*/ && (histClass.EqualTo(histClassDenom) || !fHistoList.FindObject(histClassDenom.Data())) ) continue; /// TODO: why '!optTask' was needed? else if(!histClassDenom.IsNull()) Info("DrawSame"," Denom histClass '%s' found ",histClassDenom.Data()); /// set first histogram if (i==0 && !hFirst) hFirst=h; /// print normalisation option if(optRbn||optRbnStat) Info("DrawSame"," Rebin by %d, to <%.1f%% stat. uncertainty per bin",(optRbn?rbn:0),(optRbnStat?stat*100:0)); if(optNormY||optNorm||optEvt) Info("DrawSame"," Normalize in y-axis,2D's only(%d), by int.(%d), by #events(%d)",optNormY,optNorm,optEvt); if(optSclMax) Info("DrawSame"," Scale to maximum(%d)",optSclMax); if(optCum) Info("DrawSame"," Cumulate sum of bins along x-axis, 1D's left-to-right(%d), right-to-left(%d)",!optCumR,optCumR); /// rebin, normalize, cumulate and scale spectra according to options 'rebinX','rebinStat','norm','normY','cum','cumR','events','sclMax' h->Sumw2(); if(optRbn && h->InheritsFrom( TH2::Class()) ) h=((TH2*)h)->RebinX(rbn,h->GetName()); else if( optRbn ) h->Rebin(rbn); if(optNormY && h->GetDimension()==2 && !(h->GetSumOfWeights()==0) && !optCum) PairAnalysisHelper::NormalizeSlicesY((TH2*)h); if(optRbnStat && h->GetDimension()==1) { /// rebin until stat. uncertainty is lower than 'stat' limits = PairAnalysisHelper::MakeStatBinLimits(h,stat); h=h->Rebin(limits->GetSize()-1,h->GetName(),limits->GetArray()); h->Scale(1.,"width"); // delete limits; } if(optCum) PairAnalysisHelper::Cumulate(h,optCumR,optNorm); if(optSmooth) h->Smooth(smth); /// get default histogram titles TString ytitle = h->GetYaxis()->GetTitle(); TString ztitle = h->GetZaxis()->GetTitle(); if(ytitle.Contains("{evt}")) optEvt=kFALSE; if(optNorm && !(h->GetSumOfWeights()==0) && !optCum) h=h->DrawNormalized(i>0?(optString+"same").Data():optString.Data()); if(optEvt) h->Scale(1./events); if(optSclMax) h->Scale(1./h->GetBinContent(h->GetMaximumBin())); /// modify titles according to options switch(h->GetDimension()) { case 1: if(optEvt) h->SetYTitle( (ytitle+"/N_{evt}").Data() ); if(optNorm) h->SetYTitle( (ytitle.Append(" (normalized)")).Data() ); if(optCum) h->SetYTitle( (ytitle.Prepend("cumulated ")).Data() ); if(optRatio) h->SetYTitle( "ratio" ); if(optDiv) h->SetYTitle( "ratio" ); if(optEff) h->SetYTitle( "acceptance #times efficiency" ); if(optSclMax) h->SetYTitle( (ytitle+"/N_{max}").Data() ); break; case 2: if(optEvt) h->SetZTitle( (ztitle+"/N_{evt}").Data() ); // if(optNormY) h->SetYTitle( (ytitle.Append(" (normalized)")).Data() ); if(optNorm) h->SetZTitle( (ztitle.Prepend("normalized ")).Data() ); if(optRatio) h->SetZTitle( "ratio" ); if(optDiv) h->SetZTitle( "ratio" ); if(optEff) h->SetZTitle( "acceptance #times efficiency" ); if(optSclMax) h->SetZTitle( (ztitle+"/N_{max}").Data() ); break; } /// Calculate ratios and Draw histograms if( (optEff || optRatio) && !optNorm && !optEvt && !optTask) { Info("DrawSame"," Calculate '%s' w/o normalisation and within the same task",(optEff?"efficiency":"ratio")); // TString clMC = histClass+"_MCtruth"; THashList *clDenom = (THashList*)fHistoList.FindObject( histClassDenom.Data() ); TH1 *hMC = (TH1*) h->Clone(); // needed to preserve the labeling of non-mc histogram TString histdenomMC = UserHistogram(histClassDenom.Data(),hMC); TH1 *hdenom = (TH1*) clDenom->FindObject( histdenomMC.Data() ); if(!hdenom) { Error("DrawSame","Denominator object not found"); continue; } Info("DrawSame"," Divide %s(#=%.3e) by %s(#=%.3e)",h->GetName(),h->GetEntries(),hdenom->GetName(),hdenom->GetEntries()); delete hMC; //delete the surplus object // normalize and rebin only once hdenom->Sumw2(); //why is it crashing here if(optRbnStat && (optEff || !(i%10)) ) { hdenom=hdenom->Rebin(limits->GetSize()-1,hdenom->GetName(),limits->GetArray()); hdenom->Scale(1.,"width"); } if(optRbn && (optEff || !(i%10)) ) hdenom->RebinX(rbn); if(optEvt && (optEff || !(i%10)) ) hdenom->Scale(1./events); if(!hdenom || !h->Divide(hdenom)) { Warning("DrawSame(eff/ratio)","Division failed!!!!"); continue; } } else if( optTask && (optDiv || optEff || optRatio)) { Info("DrawSame"," Calculate '%s' using different tasks",(optEff?"efficiency":(optRatio?"ratio":"divison"))); // denominators TH1* hdenom=0x0; TH1* htden =0x0; if(optEff || optRatio) { THashList *clDenom = (THashList*)fHistoList.FindObject( histClassDenom.Data() ); TH1 *hMC = (TH1*) h->Clone(); // needed to preserve the labeling of non-mc histogram TString histdenomMC = UserHistogram(histClassDenom.Data(),hMC); //delete the surplus object delete hMC; if(clDenom) { hdenom = (TH1*) clDenom->FindObject( histdenomMC.Data() ); } if(listDenom) { THashList *clTaskDen = (THashList*)listDenom->FindObject( histClassDenom.Data() ); if(clTaskDen) { // htden=(TH1*)clTaskDen->FindObject(hdenom->GetName()); htden=(TH1*)clTaskDen->FindObject( histdenomMC.Data() ); Info("DrawSame","calculate eff/ratio using task-denom: '%s' class-denom: '%s' hist-denom: '%s'", listDenom->GetName(), histClassDenom.Data(), histdenomMC.Data()); // keep only one of them, otherwise you might divide the same objects twice if(htden) hdenom=0x0; } } } // optEff || optRatio // task ratio TH1 *htnom=0x0; if(optDiv && !optEff) { Info("DrawSame"," Search for '%s' in task '%s'",histClass.Data(),listDenom->GetName()); THashList *clTaskNom = (THashList*)listDenom->FindObject( histClass.Data() ); if(!clTaskNom) continue; htnom=(TH1*)clTaskNom->FindObject(histName.Data()); if(!htnom) continue; } if(hdenom) hdenom->Sumw2(); if(htden) htden->Sumw2(); if(htnom) htnom->Sumw2(); /// check consistent bins and rebin if needed Int_t nbinsh = (h? h->GetNbinsX(): -1); Int_t nbinsd = (hdenom? hdenom->GetNbinsX(): -2); Int_t nbinstd = (htden? htden->GetNbinsX(): -3); Int_t nbinstn = (htnom? htnom->GetNbinsX(): -4); // normalize and rebin only once if(optRbn) { if(hdenom && nbinsd >nbinsh) hdenom->RebinX(rbn); if(htden && nbinstd>nbinsh) htden->RebinX(rbn); if(htnom && nbinstn>nbinsh) htnom->RebinX(rbn); } if(optRbnStat) { if(hdenom && nbinsd >nbinsh) { hdenom=hdenom->Rebin(limits->GetSize()-1,hdenom->GetName(),limits->GetArray()); hdenom->Scale(1.,"width"); } if(htden && nbinstd>nbinsh) { htden=htden->Rebin(limits->GetSize()-1,htden->GetName(),limits->GetArray()); htden->Scale(1.,"width"); } if(htnom && nbinstn>nbinsh) { htnom=htnom->Rebin(limits->GetSize()-1,htnom->GetName(),limits->GetArray()); htnom->Scale(1.,"width"); } } if(optEvt && !i ) { if(hdenom) hdenom->Scale(1./events); if(htden) htden->Scale(1./events); if(htnom) htnom->Scale(1./events); } Info("DrawSame"," Use nominator (%p,#=%.3e) and denominator 'same task & different class'(%p,#=%.3e), 'certain task & different class'(%p,#=%.3e), 'certain task & same class'(%p,#=%.3e)", h,h->GetEntries(), hdenom,(hdenom?hdenom->GetEntries():0), htden,(htden?htden->GetEntries():0), htnom,(htnom?htnom->GetEntries():0) ); // Printf("h %p (bins%d) \t hdenom %p (bins%d) \t htdenom %p (bins%d) \t htnom %p (bins%d)", // h,h->GetNbinsX(),hdenom,(hdenom?hdenom->GetNbinsX():0), // htden,(htden?htden->GetNbinsX():0),htnom,(htnom?htnom->GetNbinsX():0)); // standard ratio if(hdenom && !h->Divide(hdenom)) { Warning("DrawSame(eff/ratio)","h & denom division failed!!!!"); continue; } else if(htden && htnom && !i && !htnom->Divide(htden)) { Warning("DrawSame(eff/ratio)","task-nom/task-denom division failed!!!!"); continue; } else if(optDiv && htnom &&!h->Divide(htnom)) { Warning("DrawSame(eff/ratio)","h & task-nom division failed!!!!"); continue; } else if(htden && !h->Divide(htden)) { Warning("DrawSame(eff/ratio)","h & task-denom division failed!!!!"); continue; } } /// scale by 1./bin width, helpfull for arbitrary bin sizes if(optWdth && !optRbnStat) h->Scale(1.,"width"); /// flip bin contents of histograms to 1./content if(optOneOver){ Info("DrawSame"," Scale by 1/content"); TH1 *hOne = (TH1*) h->Clone("one"); hOne->Reset("ICSE"); for(Int_t ib=0;ib<(h->GetNbinsX()+2)*(h->GetNbinsY()+2)*(h->GetNbinsZ()+2);ib++) hOne->SetBinContent(ib,1.); if(hOne->Divide(h)) h=hOne; } /// set special geant labels if(optGeant) { Info("DrawSame"," Set GEANT bin labels"); PairAnalysisHelper::SetGEANTBinLabels(h); } /// set special pdg labels if(optPdg || optPdgClean) { Info("DrawSame"," Set PDG bin labels"); PairAnalysisHelper::SetPDGBinLabels(h,optPdgClean); } /// style histograms if not done before /// take into account reset style option 'RstSty' if(h->GetLineColor()==kBlack && !optString.Contains("col")) { // avoid color updates h->UseCurrentStyle(); PairAnalysisStyler::Style(h, i-(optRstSty?nobj:0) ); } /// some default styles for certain options if(optString.Contains("scat")) h->SetMarkerStyle(kDot); if(optString.Contains("e")) h->SetLineStyle(kSolid); // if(optString.Contains("text")) { h->SetLineColor(1); h->SetMarkerColor(1); } /// set histogram title to current histClass // h->SetName(histClass.Data()); h->SetTitle(histClass.Data()); /// add histograms to returned array if option 'goff' is active, otherwise /// draw the histogram if(optGoff) { if(optTask) h->SetTitle(Form("%s %s",GetName(),h->GetTitle())); arr->Add(h); } else if(optStack) { if(!hs) hs = new THStack("hs",Form(";%s;%s",h->GetXaxis()->GetTitle(),h->GetYaxis()->GetTitle())); hs->Add(h); } else { optString.ReplaceAll(" ",""); Info("DrawSame"," Draw object with options: '%s'",optString.Data()); // h->Draw(i>0?(optString+"same").Data():optString.Data()); if(!optSlicesY) { h->Draw(h!=hFirst?(optString+"same").Data():optString.Data()); } else if(h->GetDimension()==2) { // loop over all projections for (Int_t bin=1;binGetNbinsX();bin++) { TH1 *hp= ((TH2*)h)->ProjectionY(Form("%s_%d",h->GetName(),bin),bin,bin,"e"); if(!hp) continue; Long64_t nentries = Long64_t(hp->GetEntries()); if(nentries == 0) { delete hp; continue; } PairAnalysisStyler::Style(hp, i+(bin-1)-(optRstSty?nobj:0) ); if( optRbn ) hp->Rebin(rbn); if(optRbnStat) { /// rebin until stat. uncertainty is lower than 'stat' limits = PairAnalysisHelper::MakeStatBinLimits(hp,stat); hp=hp->Rebin(limits->GetSize()-1,hp->GetName(),limits->GetArray()); hp->Scale(1.,"width"); } if(optCum) PairAnalysisHelper::Cumulate(hp,optCumR,optNorm); if(optSmooth) hp->Smooth(smth); hp->Draw(hFirst?(optString+"same").Data():optString.Data()); } } } /// protection e.g. normalization not possible TProfile if(h && h->GetEntries()>0.) { TString ratioName=histClassDenom; TString divName=(listDenom?listDenom->GetName():""); if(optEff) divName+=(listDenom?"(MC)":""); // adapt legend name // remove reserved words TObjArray *reservedWords = fReservedWords->Tokenize(":;"); for(Int_t ir=0; irGetEntriesFast(); ir++) { // printf("histClass %s \t search for %s \n",histClass.Data(),((TObjString*)reservedWords->At(ir))->GetString().Data()); histClass.ReplaceAll( ((TObjString*)reservedWords->At(ir))->GetString(), ""); ratioName.ReplaceAll( ((TObjString*)reservedWords->At(ir))->GetString(), ""); divName.ReplaceAll( ((TObjString*)reservedWords->At(ir))->GetString(), ""); } // change default signal names to titles for(Int_t isig=0; isigGetHeader()) || src.EqualTo(leg->GetHeader()))) rpl=""; histClass.ReplaceAll(src,rpl); ratioName.ReplaceAll(src,rpl); divName.ReplaceAll(src,rpl); } // printf("histClass %s \n",histClass.Data()); // change MCtruth to MC for(Int_t isig=0; isig0) histClass.ReplaceAll( PairAnalysis::PairClassName(iptype), ""); if(ratioName.CountChar('_')>0) ratioName.ReplaceAll( PairAnalysis::PairClassName(iptype), ""); if(divName.CountChar('_')>0) divName.ReplaceAll( PairAnalysis::PairClassName(iptype), ""); } // save Dalitz and Short underscore histClass.ReplaceAll("_{Dalitz}","#{Dalitz}"); ratioName.ReplaceAll("_{Dalitz}","#{Dalitz}"); divName.ReplaceAll("_{Dalitz}","#{Dalitz}"); histClass.ReplaceAll("_{S}","#{S}"); ratioName.ReplaceAll("_{S}","#{S}"); divName.ReplaceAll("_{S}","#{S}"); // remove delimiters histClass.ReplaceAll("_"," "); ratioName.ReplaceAll("_"," "); divName.ReplaceAll("_"," "); histClass.ReplaceAll(".",""); ratioName.ReplaceAll(".",""); divName.ReplaceAll(".",""); // get Dalitz and Short back histClass.ReplaceAll("#{Dalitz}","_{Dalitz}"); ratioName.ReplaceAll("#{Dalitz}","_{Dalitz}"); divName.ReplaceAll("#{Dalitz}","_{Dalitz}"); histClass.ReplaceAll("#{S}","_{S}"); ratioName.ReplaceAll("#{S}","_{S}"); divName.ReplaceAll("#{S}","_{S}"); // remove trailing and leading spaces histClass.Remove(TString::kBoth,' '); ratioName.Remove(TString::kBoth,' '); divName.Remove(TString::kBoth,' '); //build final ratio name if(optRatio) histClass+="/"+ratioName; // delete the surplus delete reservedWords; // modify legend option TString legOpt = optString+"L"; legOpt.ReplaceAll("hist",""); legOpt.ReplaceAll("scat",""); if(legOpt.Contains("col")) legOpt=""; legOpt.ReplaceAll("z",""); legOpt.ReplaceAll("e",""); if (optTask) histClass.Prepend(Form("%s ",GetName())); if ( (optTask && optCutStep && i) || (optStack && (i-nobj)) ) histClass.Prepend("+"); if (optDiv && !optOneOver) histClass.ReplaceAll(GetName(),Form("%s/%s",GetName(),divName.Data())); if (optDiv && optOneOver) histClass.Prepend(Form("%s/",divName.Data())); if (optDet) { for (Int_t idet=kRef; idet(idet)))) histClass=PairAnalysisHelper::GetDetName(static_cast(idet)); } } // else if(nobj) histClass=""; if(optMeanX) histClass+=Form(" #LTx#GT=%.1e",h->GetMean()); if(optRmsX) histClass+=Form(" RMS(x)=%.1e",h->GetRMS()); if(optMeanY) histClass+=Form(" #LTy#GT=%.2e",h->GetMean(2)); if(optRmsY) histClass+=Form(" RMS(y)=%.2e",h->GetRMS(2)); histClass.ReplaceAll("e+00",""); // no entry for colored plots if (optLeg && leg /*&& !legOpt.Contains("col")*/) leg->AddEntry(h,histClass.Data(),legOpt.Data()); // if (leg) leg->AddEntry(h,classTable->GetName(),(optString+"L").Data()); ++i; } else if(nobj&&leg) leg->AddEntry(hFirst,"",""); //++i; // } // cleanup if(limits) delete limits; } //while loop histclass /// draw stack histogram if(optStack) { optString.ReplaceAll(" ",""); Info("DrawSame"," Draw stacked object with options: '%s'",optString.Data()); hs->Draw(nobj>0?(optString+"same").Data():optString.Data()); } /// draw legend only once /// set legend coordinates, margins, fillstyle, fontsize if (leg) { PairAnalysisStyler::SetLegendAttributes(leg,optLegFull); if(!nobj) leg->Draw(); } /// automatic axis minimum and maximum if(gPad) { Double_t max=-1e10; Double_t min=+1e10; TListIter nextObj(gPad->GetListOfPrimitives(),kIterBackward); // TObject *obj; while ((obj = nextObj())) { if(obj->InheritsFrom(TH1::Class())) { TH1* h1 = static_cast(obj); max = TMath::Max(max,PairAnalysisHelper::GetContentMaximum(h1)); //h1->GetMaximum(); Double_t tmpmax = max*(gPad->GetLogy()?5.:1.1); if(optEff) tmpmax=1.1; h1->SetMaximum(tmpmax); Double_t objmin=PairAnalysisHelper::GetContentMinimum(h1); if( gPad->GetLogy() && objmin<0.) objmin=0.5; min=TMath::Min(min,objmin); Double_t tmpmin = min*(min<0.?1.1:0.9); // if(optEff) tmpmin=0.; h1->SetMinimum(tmpmin); // Printf("after %s max%f \t min%f \t for logy %.3e >? %.3e", // h1->GetTitle(),tmpmax,tmpmin, tmpmax/(tmpmin>0.?tmpmin:1.),TMath::Power(10.,TGaxis::GetMaxDigits())); // Printf("max%f \t min%f",h1->GetMaximum(),PairAnalysisHelper::GetContentMinimum(h1)); // max=TMath::Max(max,PairAnalysisHelper::GetContentMaximum(h1)); // min=TMath::Min(min,PairAnalysisHelper::GetContentMinimum(h1));//hobj->GetBinContent(hobj->GetMinimumBin())); // Printf("max%f \t min%f",max,min); // if(!optEff) h1->SetMaximum(max*(gPad->GetLogy()?5.:1.1)); // else h1->SetMaximum(1.1); // if(!optEff) h1->SetMinimum( min*(min<0.?1.1:0.9) ); //TODO: doesnt work, why?? Negative values? /// automatically set log option labels if(gPad->GetLogy() && (tmpmax/(tmpmin>0.?tmpmin:1.) > TMath::Power(10.,TGaxis::GetMaxDigits()) || tmpminTMath::Power(10.,+TGaxis::GetMaxDigits()) )) { // if(gPad->GetLogy() && tmpmax/(tmpmin>0.?tmpmin:1.) > TMath::Power(10.,TGaxis::GetMaxDigits())) { h1->GetYaxis()->SetMoreLogLabels(kFALSE); h1->GetYaxis()->SetNoExponent(kFALSE); } Double_t tmpXmin = h1->GetXaxis()->GetXmin(); if(gPad->GetLogx() && h1->GetXaxis()->GetXmax()/(TMath::Abs(tmpXmin)<1.e-10?1.:tmpXmin) > TMath::Power(10.,TGaxis::GetMaxDigits())) { // printf("Xaxis: max%f , min%f \t ratio %.3e >? %.3e \n",h1->GetXaxis()->GetXmax(),(TMath::Abs(tmpXmin)<1.e-10?1.:tmpXmin), // h1->GetXaxis()->GetXmax()/(TMath::Abs(tmpXmin)<1.e-10?1.:tmpXmin),TMath::Power(10.,TGaxis::GetMaxDigits())); h1->GetXaxis()->SetMoreLogLabels(kFALSE); h1->GetXaxis()->SetNoExponent(kFALSE); } } } /// draw only once the default metadata if option 'meta' is active if(!nobj && optMeta && fMetaData && !gPad->GetPrimitive("meta")) { fMetaData->DrawSame(""); } /// force legend to be drawn always on top, remove multiple versions of it /// they show up when one uses the 'task' draw option if(!optTask && 0) { if(leg) leg->DrawClone(); nextObj.Reset(); Int_t ileg = 0; while ((obj = nextObj())) { if(obj->InheritsFrom(TLegend::Class())) { if(ileg>0) delete obj; ileg++; } } } } /// styling /// NOTE: this draws a copy of the first histogram //if(gPad && !optGoff) gPad->RedrawAxis(); /// remove canvas if graphics are switched off ('goff') // if(optGoff) { c->Close(); delete c; } /// clean up if(selections) delete selections; return arr; } //_____________________________________________________________________________ void PairAnalysisHistos::SetReservedWords(const char* words) { // // set reserved words // (*fReservedWords)=words; } //_____________________________________________________________________________ void PairAnalysisHistos::StoreVariables(TObject *obj, UInt_t valType[20]) { // // // if (!obj) return; if (obj->InheritsFrom(TH1::Class())) StoreVariables(static_cast(obj), valType); else if (obj->InheritsFrom(THnBase::Class())) StoreVariables(static_cast(obj), valType); return; } //_____________________________________________________________________________ void PairAnalysisHistos::StoreVariables(TH1 *obj, UInt_t valType[20]) { // // store variables in the axis (special for TProfile3D) // Int_t dim = obj->GetDimension(); // dimension correction for profiles if(obj->IsA() == TProfile::Class() || obj->IsA() == TProfile2D::Class() || obj->IsA() == TProfile3D::Class()) { dim++; } switch( dim ) { case 4: obj->SetUniqueID(valType[3]); // Tprofile3D variable case 3: obj->GetZaxis()->SetUniqueID(valType[2]); case 2: obj->GetYaxis()->SetUniqueID(valType[1]); case 1: obj->GetXaxis()->SetUniqueID(valType[0]); // fall through is intended } return; } //_____________________________________________________________________________ void PairAnalysisHistos::StoreVariables(THnBase *obj, UInt_t valType[20]) { // // store variables in the axis // obj->Sumw2(); // check for formulas and skip the rest if needed TList *list = obj->GetListOfFunctions(); if(obj->IsA()==PairAnalysisHn::Class()) list = (static_cast(obj))->GetListOfFunctions(); if(list && list->Last()) return; Int_t dim = obj->GetNdimensions(); for(Int_t it=0; itGetAxis(it)->SetUniqueID(valType[it]); obj->GetAxis(it)->SetName(Form("%s", PairAnalysisVarManager::GetValueName(valType[it]))); obj->GetAxis(it)->SetTitle(Form("%s %s", PairAnalysisVarManager::GetValueLabel(valType[it]), PairAnalysisVarManager::GetValueUnit(valType[it]))); } return; } //_____________________________________________________________________________ void PairAnalysisHistos::FillValues(TObject *obj, const Double_t *values) { // // // if (!obj) return; if (obj->InheritsFrom(TH1::Class())) FillValues(static_cast(obj), values); else if (obj->InheritsFrom(THnBase::Class())) FillValues(static_cast(obj), values); return; } //_____________________________________________________________________________ void PairAnalysisHistos::FillValues(TH1 *obj, const Double_t *values) { // // fill values for TH1 inherted classes // Int_t dim = obj->GetDimension(); Bool_t bprf = kFALSE; // UInt_t nValues = (UInt_t) PairAnalysisVarManager::kNMaxValues; UInt_t valueTypes=obj->GetUniqueID(); if (valueTypes==(UInt_t)PairAnalysisHistos::kNoAutoFill) return; Bool_t weight = (valueTypes!=kNoWeights); // check if tprofile if(obj->IsA() == TProfile::Class() || obj->IsA() == TProfile2D::Class() || obj->IsA() == TProfile3D::Class() ) bprf=kTRUE; // TO BEAUTIFY: switch off manually weighting of profile3Ds if(obj->IsA() == TProfile3D::Class()) weight=kFALSE; // check for formulas TList *list = obj->GetListOfFunctions(); TFormula *xform = dynamic_cast(list->FindObject("xFormula")); TFormula *yform = dynamic_cast(list->FindObject("yFormula")); TFormula *zform = dynamic_cast(list->FindObject("zFormula")); TFormula *pform = dynamic_cast(list->FindObject("pFormula")); TFormula *wform = dynamic_cast(list->FindObject("wFormula")); // Bool_t bform = (xform || yform || zform /*|| wform*/); // get variables from axis unique IDs UInt_t value1=obj->GetXaxis()->GetUniqueID(); UInt_t value2=obj->GetYaxis()->GetUniqueID(); UInt_t value3=obj->GetZaxis()->GetUniqueID(); UInt_t value4=obj->GetUniqueID(); // get weighting/profile var stored in the unique ID Double_t fvals[4] = { values[value1], values[value2], values[value3], values[value4] }; // use formulas to update fill values if(xform) fvals[0] = PairAnalysisHelper::EvalFormula(xform, values); if(yform) fvals[1] = PairAnalysisHelper::EvalFormula(yform, values); if(zform) fvals[2] = PairAnalysisHelper::EvalFormula(zform, values); if(wform) fvals[3] = PairAnalysisHelper::EvalFormula(wform, values); if(pform) fvals[3] = PairAnalysisHelper::EvalFormula(pform, values); // weighting overwriting for Profile3D /* // ask for inclusive trigger map variables if(value1!=PairAnalysisVarManager::kTriggerInclONL && value1!=PairAnalysisVarManager::kTriggerInclOFF && value2!=PairAnalysisVarManager::kTriggerInclONL && value2!=PairAnalysisVarManager::kTriggerInclOFF && value3!=PairAnalysisVarManager::kTriggerInclONL && value3!=PairAnalysisVarManager::kTriggerInclOFF && value4!=PairAnalysisVarManager::kTriggerInclONL && value4!=PairAnalysisVarManager::kTriggerInclOFF ) { // no trigger map variable selected */ switch ( dim ) { case 1: if(!bprf && !weight) obj->Fill(fvals[0]); // histograms else if(!bprf && weight) obj->Fill(fvals[0], fvals[3]); // weighted histograms else if(bprf && !weight) ((TProfile*)obj)->Fill(fvals[0],fvals[1]); // profiles else ((TProfile*)obj)->Fill(fvals[0],fvals[1], fvals[3]); // weighted profiles break; case 2: if(!bprf && !weight) obj->Fill(fvals[0], fvals[1]); // histograms else if(!bprf && weight) ((TH2*)obj)->Fill(fvals[0], fvals[1], fvals[3]); // weighted histograms else if(bprf && !weight) ((TProfile2D*)obj)->Fill(fvals[0], fvals[1], fvals[2]); // profiles else ((TProfile2D*)obj)->Fill(fvals[0], fvals[1], fvals[2], fvals[3]); // weighted profiles break; case 3: if(!bprf && !weight) ((TH3*)obj)->Fill(fvals[0], fvals[1], fvals[2]); // histograms else if(!bprf && weight) ((TH3*)obj)->Fill(fvals[0], fvals[1], fvals[2], fvals[3]); // weighted histograms else if(bprf && !weight) ((TProfile3D*)obj)->Fill(fvals[0], fvals[1], fvals[2], fvals[3]); // profiles else Printf(" WARNING: weighting NOT yet possible for TProfile3Ds !"); break; } /* } else { // fill inclusive trigger map variables if(weight) return; switch ( dim ) { case 1: for(Int_t i=0; i<30; i++) { if(TESTBIT((UInt_t)fvals[0],i)) obj->Fill(i); } break; case 2: if((value1==PairAnalysisVarManager::kTriggerInclOFF && value2==PairAnalysisVarManager::kTriggerInclONL) || (value1==PairAnalysisVarManager::kTriggerInclONL && value2==PairAnalysisVarManager::kTriggerInclOFF) ) { for(Int_t i=0; i<30; i++) { if((UInt_t)fvals[0]==BIT(i)) { for(Int_t i2=0; i2<30; i2++) { if((UInt_t)fvals[1]==BIT(i2)) { obj->Fill(i, i2); } // bit fired } //loop 2 }//bit fired } // loop 1 } else if(value1==PairAnalysisVarManager::kTriggerInclONL || value1==PairAnalysisVarManager::kTriggerInclOFF) { for(Int_t i=0; i<30; i++) { if(TESTBIT((UInt_t)fvals[0],i)) obj->Fill(i, fvals[1]); } } else if(value2==PairAnalysisVarManager::kTriggerInclONL || value2==PairAnalysisVarManager::kTriggerInclOFF) { for(Int_t i=0; i<30; i++) { if(TESTBIT((UInt_t)fvals[1],i)) obj->Fill(fvals[0], i); } } else //makes no sense return; break; default: return; } } //end: trigger filling */ return; } //_____________________________________________________________________________ void PairAnalysisHistos::FillValues(THnBase *obj, const Double_t *values) { // // fill values for THn inherted classes // // skip if manual filling UInt_t value4 =obj->GetUniqueID(); // weighting variable if any if (value4==(UInt_t)PairAnalysisHistos::kNoAutoFill) return; // check for formulas and skip the rest if needed TList *list = obj->GetListOfFunctions(); if(obj->IsA()==PairAnalysisHn::Class()) list = (static_cast(obj))->GetListOfFunctions(); Bool_t useFormulas = (list && list->Last()); // do weighting Bool_t weight =(value4!=kNoWeights); // fill array const Int_t dim = obj->GetNdimensions(); Double_t fill[dim]; // loop over all axes for(Int_t it=0; it(list->FindObject(formName)); fill[it] = PairAnalysisHelper::EvalFormula(form, values); } else { fill[it] = values[obj->GetAxis(it)->GetUniqueID()]; } } // fill object if(!weight) obj->Fill(fill); else obj->Fill(fill, values[value4]); return; } //_____________________________________________________________________________ void PairAnalysisHistos::FillVarArray(TObject *obj, UInt_t *valType) { // // extract variables stored in the axis (special for TProfile3D) // if (!obj) return; if (obj->InheritsFrom(TH1::Class())) { valType[0]=((TH1*)obj)->GetXaxis()->GetUniqueID(); valType[1]=((TH1*)obj)->GetYaxis()->GetUniqueID(); valType[2]=((TH1*)obj)->GetZaxis()->GetUniqueID(); valType[3]=((TH1*)obj)->GetUniqueID(); // weighting(profile) var stored in unique ID } else if (obj->InheritsFrom(THnBase::Class())) { for(Int_t it=0; it<((THn*)obj)->GetNdimensions(); it++) valType[it]=((THn*)obj)->GetAxis(it)->GetUniqueID(); } valType[19]=obj->GetUniqueID(); //weights return; } //_____________________________________________________________________________ void PairAnalysisHistos::AdaptNameTitle(TH1 *hist, const char* histClass) { // // adapt name and title of the histogram // Int_t dim = hist->GetDimension(); TString currentName = "";//hist->GetName(); TString currentTitle = "";//hist->GetTitle(); TString hclass = histClass; //get reserved class TObjArray *arr=hclass.Tokenize("_."); arr->SetOwner(); hclass=((TObjString*)arr->At(0))->GetString(); delete arr; hclass.ToLower(); Bool_t bname = (currentName.IsNull()); Bool_t btitle = (currentTitle.IsNull()); Bool_t bprf = kFALSE; if(hist->IsA() == TProfile::Class() || hist->IsA() == TProfile2D::Class() || hist->IsA() == TProfile3D::Class()) bprf=kTRUE; // tprofile options Double_t pmin=0., pmax=0.; TString option = "", calcrange=""; Bool_t bStdOpt=kTRUE; if(bprf) { switch( dim ) { case 3: option = ((TProfile3D*)hist)->GetErrorOption(); pmin = ((TProfile3D*)hist)->GetTmin(); pmax = ((TProfile3D*)hist)->GetTmax(); break; case 2: option = ((TProfile2D*)hist)->GetErrorOption(); pmin = ((TProfile2D*)hist)->GetZmin(); pmax = ((TProfile2D*)hist)->GetZmax(); break; case 1: option = ((TProfile*)hist)->GetErrorOption(); pmin = ((TProfile*)hist)->GetYmin(); pmax = ((TProfile*)hist)->GetYmax(); break; } if(option.Contains("s",TString::kIgnoreCase)) bStdOpt=kFALSE; if(pmin!=pmax) calcrange=Form("#cbar_{%+.*f}^{%+.*f}", PairAnalysisHelper::GetPrecision(pmin),pmin, PairAnalysisHelper::GetPrecision(pmax),pmax); } UInt_t varx = hist->GetXaxis()->GetUniqueID(); UInt_t vary = hist->GetYaxis()->GetUniqueID(); UInt_t varz = hist->GetZaxis()->GetUniqueID(); UInt_t varp = hist->GetUniqueID(); Bool_t weight = (varp!=kNoWeights); if(bprf && dim==3) weight=kFALSE; // no weighting for profile3D // store titles in the axis if(btitle) { TString tit = ""; /////// set NAMES hist->GetXaxis()->SetName(PairAnalysisVarManager::GetValueName(varx)); hist->GetYaxis()->SetName(PairAnalysisVarManager::GetValueName(vary)); hist->GetZaxis()->SetName(PairAnalysisVarManager::GetValueName(varz)); // adapt according to formula TFormula *xform = dynamic_cast(hist->GetListOfFunctions()->FindObject("xFormula")); TFormula *yform = dynamic_cast(hist->GetListOfFunctions()->FindObject("yFormula")); TFormula *zform = dynamic_cast(hist->GetListOfFunctions()->FindObject("zFormula")); TFormula *wform = dynamic_cast(hist->GetListOfFunctions()->FindObject("wFormula")); if(xform) { hist->GetXaxis()->SetName(PairAnalysisHelper::GetFormulaName(xform).Data()); } if(yform) { hist->GetYaxis()->SetName(PairAnalysisHelper::GetFormulaName(yform).Data()); } if(zform) { hist->GetZaxis()->SetName(PairAnalysisHelper::GetFormulaName(zform).Data()); } /////// set TITLE hist->GetXaxis()->SetTitle(PairAnalysisVarManager::GetValueLabel(varx)); hist->GetYaxis()->SetTitle(PairAnalysisVarManager::GetValueLabel(vary)); hist->GetZaxis()->SetTitle(PairAnalysisVarManager::GetValueLabel(varz)); // adapt according to formula if(xform){hist->GetXaxis()->SetTitle(PairAnalysisHelper::GetFormulaTitle(xform).Data()); } if(yform){hist->GetYaxis()->SetTitle(PairAnalysisHelper::GetFormulaTitle(yform).Data()); } if(zform){hist->GetZaxis()->SetTitle(PairAnalysisHelper::GetFormulaTitle(zform).Data()); } // profile axis if(bprf && dim<3) { TAxis *ax = 0x0; switch( dim ) { case 2: ax = hist->GetZaxis(); break; case 1: ax = hist->GetYaxis(); break; } tit=ax->GetTitle(); tit.Prepend( (bStdOpt ? "#LT" : "RMS(") ); tit.Append ( (bStdOpt ? "#GT" : ")") ); // TODO: activate --> tit.Append ( Form("_{%ss}",hclass.Data())); tit.Append ( calcrange.Data() ); ax->SetTitle(tit.Data()); } // append the units for all axes (except formula) tit=Form("%s %s",hist->GetXaxis()->GetTitle(), PairAnalysisVarManager::GetValueUnit(varx)); if(!xform) hist->GetXaxis()->SetTitle(tit.Data()); tit=Form("%s %s",hist->GetYaxis()->GetTitle(), PairAnalysisVarManager::GetValueUnit(vary)); if(!yform) hist->GetYaxis()->SetTitle(tit.Data()); tit=Form("%s %s",hist->GetZaxis()->GetTitle(), PairAnalysisVarManager::GetValueUnit(varz)); if(!zform) hist->GetZaxis()->SetTitle(tit.Data()); // overwrite titles with hist class if needed if(!bprf) { switch( dim ) { case 1: hist->GetYaxis()->SetTitle(Form("%ss",hclass.Data())); break; case 2: hist->GetZaxis()->SetTitle(Form("%ss",hclass.Data())); break; } } // weighted axis (maximal 2 dimensional) if(weight) { TAxis *ax = hist->GetYaxis(); if(dim==2)ax = hist->GetZaxis(); tit=PairAnalysisVarManager::GetValueLabel(varp); if(wform) { tit=PairAnalysisHelper::GetFormulaTitle(wform); } ax->SetTitle(Form("%s weighted %s", tit.Data(), ax->GetTitle())); } // create an unique name TFormula *pform = dynamic_cast(hist->GetListOfFunctions()->FindObject("pFormula")); if(bname) switch(dim) { case 3: currentName+=Form("%s_",hist->GetXaxis()->GetName()); currentName+=Form("%s_",hist->GetYaxis()->GetName()); currentName+=Form("%s", hist->GetZaxis()->GetName()); if(bprf && !pform) currentName+=Form("-%s%s",PairAnalysisVarManager::GetValueName(varp),(bStdOpt ? "avg" : "rms")); else if (bprf) currentName+=Form("-%s%s",PairAnalysisHelper::GetFormulaName(pform).Data(),(bStdOpt ? "avg" : "rms")); if(weight&&!bprf) currentName+=Form("-wght%s",PairAnalysisVarManager::GetValueName(varp)); break; case 2: currentName+=Form("%s_",hist->GetXaxis()->GetName()); currentName+=Form("%s", hist->GetYaxis()->GetName()); if(bprf) currentName+=Form("-%s%s",hist->GetZaxis()->GetName(),(bStdOpt ? "avg" : "rms")); if(weight && !wform) currentName+=Form("-wght%s",PairAnalysisVarManager::GetValueName(varp)); else if(weight && wform) currentName+=Form("-wght%s",PairAnalysisHelper::GetFormulaName(wform).Data()); break; case 1: currentName+=Form("%s",hist->GetXaxis()->GetName()); if(bprf) currentName+=Form("-%s%s",hist->GetYaxis()->GetName(),(bStdOpt ? "avg" : "rms")); if(weight && !wform) currentName+=Form("-wght%s",PairAnalysisVarManager::GetValueName(varp)); else if(weight && wform) currentName+=Form("-wght%s",PairAnalysisHelper::GetFormulaName(wform).Data()); break; } // to differentiate btw. leg and pair histos // if(!strcmp(histClass,"Pair")) currentName.Prepend("p"); if(hclass.Contains("pair")) currentName.Prepend("p"); hist->SetName(currentName.Data()); } } //_____________________________________________________________________________ void PairAnalysisHistos::AdaptNameTitle(THnBase *hist, const char* histClass) { // // adapt name and title of the histogram // Int_t dim = hist->GetNdimensions(); TString currentName = "";//hist->GetName(); TString currentTitle = "";//hist->GetTitle(); TString hclass = histClass; //get reserved class TObjArray *arr=hclass.Tokenize("_."); arr->SetOwner(); hclass=((TObjString*)arr->At(0))->GetString(); delete arr; hclass.ToLower(); Bool_t bname = (currentName.IsNull()); Bool_t btitle = (currentTitle.IsNull()); TList *list = (static_cast(hist))->GetListOfFunctions(); // store titles in the axis if(btitle) { TString tit = ""; // adapt according to formula for(Int_t it=0; it(list->FindObject(formName)); hist->GetAxis(it)->SetName(PairAnalysisHelper::GetFormulaName(form).Data()); // adapt according to formula hist->GetAxis(it)->SetTitle(PairAnalysisHelper::GetFormulaTitle(form).Data()); } } // create an unique name if(bname) { currentName=hist->GetName(); if(hclass.Contains("pair")) currentName.Prepend("p"); hist->SetName(currentName.Data()); } }