PWGHF/vertexingHF/macros/DrawQAoutput.C

   1 #if !defined(__CINT__) || defined(__MAKECINT__)
   2 #include <Riostream.h>
   3 #include <fstream>
   4 #include <TFile.h>
   5 #include <TString.h>
   6 #include <TH2F.h>
   7 #include <TH1F.h>
   8 #include <TF1.h>
   9 #include <TGraph.h>
  10 #include <TDirectoryFile.h>
  11 #include <TList.h>
  12 #include <TCanvas.h>
  13 #include <TLegend.h>
  14 #include <TPaveText.h>
  15 #include <TPaveStats.h>
  16 #include <TStyle.h>
  17 #include <TClass.h>
  18 #include <TDatabasePDG.h>
  19 #include <TParameter.h>
  20 #include <AliCounterCollection.h>
  21 #include <AliRDHFCuts.h>
  22 #endif
  23
  24 /* $Id$ */
  25
  26 TString *periodsname;
  27
  28 //read the file and take list and stat
  29
  30 Bool_t ReadFile(TList* &list,TH1F* &hstat, TString listname,TString partname,TString path="./",TString filename=/*"AnalysisResults.root"*/"PWG3histograms.root", TString dirname="PWG3_D2H_QA");
  31 Bool_t ReadFileMore(TList* &list,TH1F* &hstat, AliRDHFCuts* &cutobj, TString listname,TString partname,TString path="./",TString filename=/*"AnalysisResults.root"*/"PWG3histograms.root", TString dirname="PWG3_D2H_QA");
  32 void SuperimposeBBToTPCSignal(Int_t period /*0=LHC10bc, 1=LHC10d, 2=LHC10h*/,TCanvas* cpid, Int_t set);
  33 void TPCBetheBloch(Int_t set);
  34 Bool_t ReadFilesForCompilation(TString inputtextfile, TList**& lists, Int_t&numb, TString*& legend);
  35 void CompilationEvSelection(Int_t n,TList** lists, TString* legend);
  36 void CompilationTrackSelection(Int_t n,TList** lists, TString* legend);
  37
  38 Bool_t ReadFile(TList* &list,TH1F* &hstat, TString listname,TString partname,TString path,TString filename, TString dirname){
  39
  40   TString hstatname="nEntriesQA", cutobjname="";
  41   filename.Prepend(path);
  42   listname+=partname;
  43   hstatname+=partname;
  44
  45   TFile* f=new TFile(filename.Data());
  46   if(!f->IsOpen()){
  47     cout<<filename.Data()<<" not found"<<endl;
  48     return kFALSE;
  49   }
  50   TDirectoryFile* dir=(TDirectoryFile*)f->Get(dirname);
  51   if(!dir){
  52     cout<<dirname.Data()<<" not found in "<<filename.Data()<<endl;
  53     f->ls();
  54     return kFALSE;
  55   }
  56
  57   list=(TList*)dir->Get(listname);
  58   if(!list){
  59     cout<<"List "<<listname.Data()<<" not found"<<endl;
  60     dir->ls();
  61     return kFALSE;
  62   }
  63
  64   hstat=(TH1F*)dir->Get(hstatname);
  65   if(!hstat){
  66     cout<<hstatname.Data()<<" not found"<<endl;
  67     return kFALSE;
  68   }
  69
  70   return kTRUE;
  71 }
  72
  73 Bool_t ReadFileMore(TList* &list,TH1F* &hstat, AliRDHFCuts* &cutobj, TString listname,TString partname,TString path,TString filename,TString dirname){
  74
  75   TString hstatname="nEntriesQA", cutobjname="";
  76   filename.Prepend(path);
  77   listname+=partname;
  78   hstatname+=partname;
  79
  80   if(partname.Contains("Dplus")) cutobjname="AnalysisCuts";//"DplustoKpipiCutsStandard";
  81   else{
  82     if(partname.Contains("D0")) cutobjname="D0toKpiCutsStandard";//"D0toKpiCuts"
  83     else{
  84       if(partname.Contains("Dstar")) cutobjname="DStartoKpipiCuts";
  85       else{
  86         if(partname.Contains("Ds")) cutobjname="DstoKKpiCuts";
  87         else{
  88           if(partname.Contains("D04")) cutobjname="D0toKpipipiCuts";
  89           else{
  90             if(partname.Contains("Lc")) cutobjname="LctopKpiAnalysisCuts";
  91           }
  92         }
  93       }
  94     }
  95   }
  96
  97   TFile* f=new TFile(filename.Data());
  98   if(!f->IsOpen()){
  99     cout<<filename.Data()<<" not found"<<endl;
 100     return kFALSE;
 101   }
 102   TDirectoryFile* dir=(TDirectoryFile*)f->Get(dirname);
 103   if(!dir){
 104     cout<<dirname.Data()<<" not found  in "<<filename.Data()<<endl;
 105     return kFALSE;
 106   }
 107
 108   list=(TList*)dir->Get(listname);
 109   if(!list){
 110     cout<<"List "<<listname.Data()<<" not found"<<endl;
 111     dir->ls();
 112     return kFALSE;
 113   }
 114
 115   hstat=(TH1F*)dir->Get(hstatname);
 116   if(!hstat){
 117     cout<<hstatname.Data()<<" not found"<<endl;
 118     return kFALSE;
 119   }
 120
 121   cutobj=(AliRDHFCuts*)dir->Get(cutobjname);
 122   if(!cutobj){
 123     cout<<cutobjname.Data()<<" not found"<<endl;
 124     return kFALSE;
 125   }
 126
 127   return kTRUE;
 128 }
 129
 130 //draw "track related" histograms (list "outputTrack")
 131 void DrawOutputTrack(TString partname="D0",TString textleg="",TString path="./", Bool_t superimpose=kFALSE, TString suffixdir="",TString filename=/*"AnalysisResults.root"*/"PWG3histograms.root", Bool_t normint=kTRUE /*normalize at integral, or at nevents*/){
 132   gStyle->SetCanvasColor(0);
 133   gStyle->SetTitleFillColor(0);
 134   gStyle->SetStatColor(0);
 135   gStyle->SetPalette(1);
 136
 137   TString listname="outputTrack",name1="",name2="",path2="",filename2="PWG3histograms.root",legtext1="",legtext2="",suffixdir2="";
 138   TString tmp="y";
 139
 140   if(superimpose){
 141     cout<<"Enter the names:\n>";
 142     cin>>name1;
 143     cout<<">";
 144     cin>>name2;
 145     cout<<"Are they in the same output file? (y/n)"<<endl;
 146     cin>>tmp;
 147     if(tmp=="n"){
 148       cout<<"Path: \n";
 149       cout<<">";
 150       cin>>path2;
 151       cout<<"Filename: "<<endl;
 152       cout<<">";
 153       cin>>filename2;
 154       cout<<"Suffix in the dir name, if any (otherwhise write no)\n>";
 155       cin>>suffixdir2;
 156       if(suffixdir2=="no") suffixdir2="";
 157       cout<<"Text in the legend:\n 1)";
 158       cin>>legtext1;
 159       cout<<"2)";
 160       cin>>legtext2;
 161     }
 162
 163   }
 164
 165   Float_t nevents,nevents2;
 166   TList* list;
 167   TH1F * hstat;
 168   TString dirname="PWG3_D2H_QA";
 169   //dirname+=suffixdir;
 170   Bool_t isRead=ReadFile(list,hstat,listname,Form("%s%s",partname.Data(),name1.Data()),path,filename,Form("%s%s",dirname.Data(),suffixdir.Data()));
 171   if(!isRead) return;
 172   if(!list || !hstat){
 173     cout<<":-( null pointers..."<<endl;
 174     return;
 175   }
 176   nevents=hstat->GetBinContent(1);
 177   TPaveText *pvtxt=new TPaveText(0.6,0.6,0.9,0.9,"NDC");
 178   pvtxt->SetBorderSize(0);
 179   pvtxt->SetFillStyle(0);
 180   pvtxt->AddText(legtext1);
 181
 182   TList* llist;
 183   TH1F* hhstat;
 184   if(superimpose){
 185     isRead=ReadFile(llist,hhstat,listname,Form("%s%s",partname.Data(),name2.Data()),path2,filename2,Form("%s%s",dirname.Data(),suffixdir2.Data()));
 186     if(!isRead) return;
 187     if(!llist || !hhstat){
 188       cout<<":-( null pointers..."<<endl;
 189       return;
 190     }
 191     nevents2=hhstat->GetBinContent(1);
 192     TText *redtext=pvtxt->AddText(legtext2);
 193     redtext->SetTextColor(kRed);
 194     hhstat->Scale(hstat->Integral()/hhstat->Integral());
 195
 196   }
 197
 198   for(Int_t i=0;i<list->GetEntries();i++){
 199     TH1F* h=(TH1F*)list->At(i);
 200     TH1F* hh=0x0;
 201     TH1F* hr=0x0;
 202     if(superimpose){
 203       hh=(TH1F*)llist->At(i);
 204       hr=(TH1F*)hh->Clone(Form("%s_ratio",hh->GetName()));
 205       if(hh && TMath::Abs(hh->Integral()) > 1e-6)hh->Scale(h->Integral()/hh->Integral());
 206     }
 207     if(!h || (superimpose && !hh)){
 208       cout<<"Histogram "<<i<<" not found"<<endl;
 209       continue;
 210     }
 211     if(superimpose){
 212       if(normint) hh->Scale(h->Integral()/hh->Integral());
 213       else hh->Scale(nevents/nevents2);
 214       hhstat->SetLineColor(kRed);
 215       hh->SetLineColor(kRed);
 216       hr->Divide(h);
 217     }
 218
 219     TCanvas* c=new TCanvas(Form("c%s",h->GetName()),h->GetName());
 220     c->cd();
 221     c->SetGrid();
 222     TString hname=h->GetName();
 223     if(!hname.Contains("nCls")){
 224       c->SetLogy();
 225       if(hname.Contains("Layer")){
 226         for(Int_t ibin=1;ibin<=h->GetNbinsX();ibin++){
 227           h->GetXaxis()->SetBinLabel(ibin+1,Form("%d",ibin));
 228         }
 229         h->GetXaxis()->SetLabelSize(0.06);
 230         h->GetXaxis()->SetRangeUser(0,6); //comment to see ntracks!
 231       }
 232       //h->SetMinimum(1);
 233       h->Draw();
 234       if(superimpose)
 235         {
 236           hh->Draw("sames");
 237           TCanvas* c2=new TCanvas(Form("c2%s",h->GetName()),h->GetName());
 238           c2->cd();
 239           c2->SetGrid();
 240           hr->Draw();
 241           c2->SaveAs(Form("%s%s%s%sRatio.png",c->GetName(),name1.Data(),name2.Data(),textleg.Data()));
 242
 243         }
 244     } else {
 245       h->Draw("htext0");
 246       if(superimpose)hh->Draw("htext0sames");
 247     }
 248     c->cd();
 249     pvtxt->Draw();
 250     c->SaveAs(Form("%s%s%s%s.png",c->GetName(),name1.Data(),name2.Data(),textleg.Data()));
 251     c->SaveAs(Form("%s%s%s%s.eps",c->GetName(),name1.Data(),name2.Data(),textleg.Data()));
 252   }
 253
 254   TCanvas* cst=new TCanvas("cst","Stat");
 255   cst->SetGridy();
 256   cst->cd();
 257   hstat->Draw("htext0");
 258   if(superimpose) {
 259     hhstat->Draw("htext0sames");
 260     pvtxt->Draw();
 261   }
 262   cst->SaveAs(Form("%s%s.png",hstat->GetName(),textleg.Data()));
 263   cst->SaveAs(Form("%s%s.eps",hstat->GetName(),textleg.Data()));
 264
 265   TH1F* hd0fb4=(TH1F*)list->FindObject("hd0TracksFilterBit4");
 266   TH1F* hd0SPD1=(TH1F*)list->FindObject("hd0TracksSPDin");
 267   TH1F* hd0SPDany=(TH1F*)list->FindObject("hd0TracksSPDany");
 268   TH1F* hd0TPCITScuts=(TH1F*)list->FindObject("hd0TracksTPCITSSPDany");
 269   TH1F* hhd0fb4=0x0; TH1F* hhd0SPD1=0x0; TH1F* hhd0SPDany=0x0; TH1F* hhd0TPCITScuts=0x0;
 270   if(superimpose){
 271     hhd0fb4=(TH1F*)llist->FindObject("hd0TracksFilterBit4");
 272     hhd0SPD1=(TH1F*)llist->FindObject("hd0TracksSPDin");
 273     hhd0SPDany=(TH1F*)llist->FindObject("hd0TracksSPDany");
 274     hhd0TPCITScuts=(TH1F*)llist->FindObject("hd0TracksTPCITSSPDany");
 275
 276   }
 277   if(hd0fb4 && hd0SPD1 && hd0SPDany && hd0TPCITScuts){
 278     TCanvas* ctrsel=new TCanvas("ctrsel","Track Sel");
 279     ctrsel->SetLogy();
 280     hd0SPD1->SetLineColor(kCyan+3);
 281     hd0SPD1->Draw();
 282     ctrsel->Update();
 283     TPaveStats *st1=(TPaveStats*)hd0SPD1->GetListOfFunctions()->FindObject("stats");
 284     st1->SetTextColor(kCyan+3);
 285     st1->SetY1NDC(0.71);
 286     st1->SetY2NDC(0.9);
 287     hd0SPDany->SetLineColor(4);
 288     hd0SPDany->Draw("sames");
 289     ctrsel->Update();
 290     TPaveStats *st2=(TPaveStats*)hd0SPDany->GetListOfFunctions()->FindObject("stats");
 291     st2->SetY1NDC(0.51);
 292     st2->SetY2NDC(0.7);
 293     st2->SetTextColor(4);
 294     hd0fb4->SetLineColor(2);
 295     hd0fb4->Draw("sames");
 296     ctrsel->Update();
 297     TPaveStats *st3=(TPaveStats*)hd0fb4->GetListOfFunctions()->FindObject("stats");
 298     st3->SetY1NDC(0.31);
 299     st3->SetY2NDC(0.5);
 300     st3->SetTextColor(2);
 301     hd0TPCITScuts->SetLineColor(kGreen+1);
 302     hd0TPCITScuts->Draw("sames");
 303     ctrsel->Update();
 304     TPaveStats *st4=(TPaveStats*)hd0TPCITScuts->GetListOfFunctions()->FindObject("stats");
 305     st4->SetY1NDC(0.71);
 306     st4->SetY2NDC(0.9);
 307     st4->SetX1NDC(0.55);
 308     st4->SetX2NDC(0.75);
 309     st4->SetTextColor(kGreen+1);
 310     ctrsel->Modified();
 311     TLegend* leg=new TLegend(0.15,0.5,0.45,0.78);
 312     leg->SetFillStyle(0);
 313     leg->SetBorderSize(0);
 314     leg->AddEntry(hd0SPD1,"kITSrefit+SPD inner","L");
 315     leg->AddEntry(hd0SPDany,"kITSrefit+SPD any","L");
 316     leg->AddEntry(hd0TPCITScuts,"TPC+ITS cuts+SPD any","L");
 317     leg->AddEntry(hd0fb4,"Filter Bit 4","L");
 318     leg->Draw();
 319
 320     if(superimpose && hhd0fb4 && hhd0SPD1 && hhd0SPDany && hhd0TPCITScuts){
 321       hhd0SPD1->SetLineColor(kCyan+3);
 322       hhd0SPD1->SetStats(0);
 323       hhd0SPD1->SetLineStyle(3);
 324       hhd0SPDany->SetLineColor(4);
 325       hhd0SPDany->SetStats(0);
 326       hhd0SPDany->SetLineStyle(3);
 327       hhd0TPCITScuts->SetLineColor(kGreen+1);
 328       hhd0TPCITScuts->SetStats(0);
 329       hhd0TPCITScuts->SetLineStyle(3);
 330       hhd0fb4->SetLineColor(2);
 331       hhd0fb4->SetStats(0);
 332       hhd0fb4->SetLineStyle(3);
 333       ctrsel->cd();
 334       hhd0SPD1->Draw("sames");
 335       hhd0SPDany->Draw("sames");
 336       hhd0TPCITScuts->Draw("sames");
 337       hhd0fb4->Draw("sames");
 338
 339     }
 340     ctrsel->SaveAs("ImpactParameterTrackSel.eps");
 341     ctrsel->SaveAs("ImpactParameterTrackSel.png");
 342
 343   }
 344 }
 345
 346 //draw "pid related" histograms (list "outputPID")
 347 //period=-999 to draw the pull instead of the cut
 348 void DrawOutputPID(TString partname="D0", Int_t mode=0/*0=with pull, 1=with nsigma*/,TString textleg="",TString path="./",TString suffixdir="", TString filename="AnalysisResults.root"){
 349   gStyle->SetCanvasColor(0);
 350   gStyle->SetTitleFillColor(0);
 351   gStyle->SetOptStat(0);
 352   gStyle->SetPalette(1);
 353
 354   Int_t period=2 ,set=0;
 355   if(mode==1){
 356     cout<<"Choose period: \n-LHC10h -> 2;\n-LHC10de -> 1;\n-LHC10bc -> 0"<<endl;
 357     cin>>period;
 358     if(period>0){
 359       cout<<"Choose set: "<<endl;
 360       if(period==2) cout<<"-pass1 -> 0;\n-pass2 -> 1"<<endl;
 361       cin>>set;
 362     }
 363   }
 364
 365   TString listname="outputPid";
 366   TString dirname="PWG3_D2H_QA";
 367   dirname+=suffixdir;
 368
 369   TList* list;
 370   TH1F * hstat;
 371   //needed only for mode 1
 372   AliRDHFCuts* cutobj;
 373   AliAODPidHF* aodpid;
 374   Double_t nsigmaTOF=0;
 375   Double_t nsigmaTPC[3]={},plimTPC[2]={};
 376
 377   if(mode==1){
 378     Bool_t isRead=ReadFileMore(list,hstat,cutobj,listname,partname+suffixdir,path,filename,dirname);
 379     if(!isRead) return;
 380     if(!list || !hstat){
 381       cout<<":-( null pointers..."<<endl;
 382       return;
 383     }
 384     aodpid=(AliAODPidHF*)cutobj->GetPidHF();
 385     if(!aodpid){
 386       cout<<"PidHF object not found! cannot get the nsigma values"<<endl;
 387       return;
 388     }
 389     nsigmaTOF=aodpid->GetSigma(3);
 390
 391     nsigmaTPC[0]=aodpid->GetSigma(0);
 392     nsigmaTPC[1]=aodpid->GetSigma(1);
 393     nsigmaTPC[2]=aodpid->GetSigma(2);
 394     aodpid->GetPLimit(plimTPC);
 395
 396   }else{
 397     Bool_t isRead=ReadFile(list,hstat,listname,partname+suffixdir,path,filename,dirname);
 398     if(!isRead) return;
 399     if(!list || !hstat){
 400       cout<<":-( null pointers..."<<endl;
 401       return;
 402     }
 403   }
 404
 405
 406   TPaveText *txtsigmaTOF=new TPaveText(0.1,0.65,0.5,0.9,"NDC");
 407   txtsigmaTOF->SetBorderSize(0);
 408   txtsigmaTOF->SetFillStyle(0);
 409   txtsigmaTOF->AddText(Form("nsigmacut from cutobj = %.1f",nsigmaTOF));
 410   TLine lTOF;
 411   lTOF.SetLineColor(kMagenta+1);
 412   lTOF.SetLineStyle(2);
 413   lTOF.SetLineWidth(3);
 414
 415   TPaveText *txtsigmaTPC=new TPaveText(0.3,0.6,0.6,0.9,"NDC");
 416   txtsigmaTPC->SetBorderSize(0);
 417   txtsigmaTPC->SetFillStyle(0);
 418   txtsigmaTPC->AddText("nsigmacut from cutobj \n");
 419   txtsigmaTPC->AddText(Form("p < %.1f : %.1f \n",plimTPC[0],nsigmaTPC[0]));
 420   txtsigmaTPC->AddText(Form("%.1f < p < %.1f : %.1f \n",plimTPC[0],plimTPC[1],nsigmaTPC[1]));
 421   txtsigmaTPC->AddText(Form("p > %.1f : %.1f \n",plimTPC[1],nsigmaTPC[2]));
 422   TLine lTPC;
 423   lTPC.SetLineColor(kMagenta+1);
 424   lTPC.SetLineStyle(2);
 425   lTPC.SetLineWidth(3);
 426
 427   // TCanvas *ctest=new TCanvas("text","Test text");
 428   // ctest->cd();
 429   // txtsigmaTPC->Draw();
 430   // txtsigmaTOF->Draw();
 431
 432
 433   for(Int_t i=0;i<list->GetEntries();i++){
 434     TClass* objtype=list->At(i)->IsA();
 435     TString tpname=objtype->GetName();
 436
 437     if(tpname=="TH1F"){
 438       TH1F* h=(TH1F*)list->At(i);
 439
 440       if(!h){
 441         cout<<"Histogram "<<i<<" not found"<<endl;
 442         continue;
 443       }
 444       //h->Scale(1./h->Integral("width"));
 445       TCanvas* c=new TCanvas(Form("c%s",h->GetName()),h->GetName());
 446       c->SetLogz();
 447       c->cd();
 448       h->Draw();
 449
 450       //write
 451       c->SaveAs(Form("%s%s.png",h->GetName(),textleg.Data()));
 452       c->SaveAs(Form("%s%s.eps",h->GetName(),textleg.Data()));
 453       TFile* fout=new TFile(Form("%s.root",h->GetName()),"recreate");
 454       fout->cd();
 455       c->Write();
 456     }
 457
 458     if(tpname=="TH2F"){
 459       TH2F* h=(TH2F*)list->At(i);
 460
 461       if(!h){
 462         cout<<"Histogram "<<i<<" not found"<<endl;
 463         continue;
 464       }
 465       TString hname=h->GetName();
 466       h->Sumw2();
 467       if(h->Integral("width")==0) {cout<<"Empty histogram, skip\n"; continue;}
 468       h->Scale(1./h->Integral("width"));
 469
 470       Double_t maxzaxis=h->GetBinContent(h->GetMaximumBin());
 471       Double_t minzaxis=h->GetBinContent(h->GetMinimumBin());
 472       printf("Minimum = %f, maximum = %f\n",minzaxis,maxzaxis);
 473       TH2F* hallzrange=(TH2F*)h->Clone(Form("%swholez",hname.Data()));
 474       hallzrange->SetAxisRange(1e-07,maxzaxis,"Z");
 475       //hallzrange->SetAxisRange(minzaxis,maxzaxis,"Z");
 476
 477       TCanvas* cwholez=new TCanvas(Form("c%swholez",hname.Data()),Form("%s down to lowest z",hname.Data()));
 478       cwholez->SetLogz();
 479       hallzrange->Draw("colz");
 480       cwholez->SaveAs(Form("%swholez.png",h->GetName()));
 481       cwholez->SaveAs(Form("%swholez.eps",h->GetName()));
 482
 483       if(hname.Contains("hTOFtimeKaonHyptime")){
 484         TCanvas* cz=new TCanvas(Form("c%szoom",hname.Data()),Form("%szoom",hname.Data()));
 485         cz->SetLogz();
 486         TH2F* hz=(TH2F*)h->Clone(Form("%sz",hname.Data()));
 487         hz->Draw("colz");
 488         hz->SetAxisRange(-1500,1500,"Y");
 489         hz->SetAxisRange(0.,5.,"X");
 490         //write
 491         cz->SaveAs(Form("%szoom.png",h->GetName()));
 492         cz->SaveAs(Form("%szoom.eps",h->GetName()));
 493       }
 494
 495       TCanvas* c=new TCanvas(Form("c%s",hname.Data()),hname.Data());
 496       c->SetLogz();
 497       //c->SetLogx();
 498       TCanvas* c2=new TCanvas(Form("c2%s",hname.Data()),hname.Data());
 499       c2->SetLogz();
 500
 501       c->cd();
 502       h->DrawClone("colz");
 503
 504       if (hname.Contains("Sig") || hname.Contains("sigma"))h->SetAxisRange(-5,5,"Y");
 505       c2->cd();
 506       //if (hname.Contains("TOFtime"))h->SetAxisRange(-1500,1500,"Y");
 507       h->SetAxisRange(0.,5.,"X");
 508
 509       h->Draw("colz");
 510
 511       //TCanvas *test=new TCanvas("test","test");
 512       if(mode==0){
 513         //mean and pull, code from Jens Wiechula
 514         TF1 fg("fg","gaus",-2.,2.); // fit range +- 2 sigma
 515         TLine l;
 516         TObjArray arr;
 517
 518         //h->Draw("colz");
 519         fg.SetParameters(1,0,1);
 520         h->FitSlicesY(&fg,0,-1,0,"NQR",&arr);
 521
 522         TH1 *hM=(TH1*)arr.At(1);
 523         hM->SetMarkerStyle(20);
 524         hM->SetMarkerSize(.5);
 525         hM->DrawClone("sames");
 526
 527         TH1 *hS=(TH1*)arr.At(2);
 528         hS->SetMarkerStyle(20);
 529         hS->SetMarkerSize(.5);
 530         hS->SetMarkerColor(kRed);
 531         hS->SetLineColor(kRed);
 532         hS->DrawClone("same");
 533
 534         l.SetLineColor(kBlack);
 535         l.DrawLine(.2,0,20,0);
 536         l.SetLineColor(kRed);
 537         l.DrawLine(.2,1,20,1);
 538
 539       }else{ //mode 1
 540
 541         if(hname.Contains("TOFsigma")) {
 542
 543           c->cd();
 544           txtsigmaTOF->Draw();
 545           lTOF.DrawLine(.2,nsigmaTOF,20,nsigmaTOF);
 546           lTOF.DrawLine(.2,-1*nsigmaTOF,4.,-1*nsigmaTOF);
 547
 548         }
 549
 550
 551         if(hname.Contains("TPCsigma")){
 552
 553           c->cd();
 554           txtsigmaTPC->Draw();
 555           lTPC.DrawLine(0.,nsigmaTPC[0],plimTPC[0],nsigmaTPC[0]);
 556           lTPC.DrawLine(plimTPC[0],nsigmaTPC[1],plimTPC[1],nsigmaTPC[1]);
 557           lTPC.DrawLine(plimTPC[1],nsigmaTPC[2],4,nsigmaTPC[2]);
 558           lTPC.DrawLine(0.,-1*nsigmaTPC[0],plimTPC[0],-1*nsigmaTPC[0]);
 559           lTPC.DrawLine(plimTPC[0],-1*nsigmaTPC[1],plimTPC[1],-1*nsigmaTPC[1]);
 560           lTPC.DrawLine(plimTPC[1],-1*nsigmaTPC[2],4,-1*nsigmaTPC[2]);
 561         }
 562
 563         if(hname.Contains("TPCsigvsp")){
 564           SuperimposeBBToTPCSignal(period,c,set);
 565         }
 566       }
 567
 568       //write
 569       c->SaveAs(Form("%s%d.png",h->GetName(),mode));
 570       c->SaveAs(Form("%s%d.eps",h->GetName(),mode));
 571       c2->SaveAs(Form("%s2%d.png",h->GetName(),mode));
 572       c2->SaveAs(Form("%s2%d.eps",h->GetName(),mode));
 573
 574       TFile* fout=new TFile(Form("%s%d.root",h->GetName(),mode),"recreate");
 575       fout->cd();
 576       c->Write();
 577       c2->Write();
 578     }
 579   }
 580 }
 581
 582 void SuperimposeBBToTPCSignal(Int_t period /*0=LHC10bc, 1=LHC10d, 2=LHC10h*/,TCanvas* cpid,Int_t set /*see below*/){
 583
 584   TFile* fBethe=new TFile("BetheBlochTPC.root");
 585   if(!fBethe->IsOpen()){
 586     TPCBetheBloch(set);
 587     fBethe=new TFile("BetheBlochTPC.root");
 588   }
 589   const Int_t npart=4;
 590   TString partnames[npart]={"Kaon","Pion","Electron","Proton"};
 591   for(Int_t ipart=0;ipart<npart;ipart++){
 592     TString grname=Form("%sP%d",partnames[ipart].Data(),period);
 593     TGraph* gr=(TGraph*)fBethe->Get(grname);
 594     cpid->cd();
 595     gr->SetLineColor(1);
 596     gr->SetLineWidth(2);
 597     gr->Draw("L");
 598   }
 599
 600   //cpid->SaveAs(Form("%sBB.png",hname.Data()));
 601 }
 602
 603 //draw and save Bethe Bloch from TPC in different periods
 604 void TPCBetheBloch(Int_t set){
 605   gStyle->SetOptTitle(0);
 606   gStyle->SetCanvasColor(0);
 607
 608   AliTPCPIDResponse *tpcResp=new AliTPCPIDResponse();
 609
 610   const Int_t npart=4;
 611   //Double_t masses[npart]={TDatabasePDG::Instance()->GetParticle(321)->Mass()/*Kaon*/,TDatabasePDG::Instance()->GetParticle(211)->Mass()/*Pion*/,TDatabasePDG::Instance()->GetParticle(11)->Mass()/*Electron*/,TDatabasePDG::Instance()->GetParticle(2212)->Mass()/*Proton*/};
 612   TString partnames[npart]={"Kaon","Pion","Electron","Proton"};
 613   //printf("%s = %.4f,%s = %.4f,%s = %.4f\n",partnames[0].Data(),masses[0],partnames[1].Data(),masses[1],partnames[2].Data(),masses[2]);
 614   TCanvas *cBethe=new TCanvas("cBethe","Bethe Bloch K pi e p");
 615   Int_t nperiods=4; //LHC10b+c, LHC10d, LHC10h, MC
 616   Double_t alephParameters[5]={};
 617   Int_t nsets=1/*LHC10bc*/+2/*LHC10de*/+2/*LHC10h*/+3/*MC*/;
 618
 619   periodsname=new TString[nsets];
 620   cout<<"Creating the file of the Bethe Bloch"<<endl;
 621   TFile* fout=new TFile("BetheBlochTPC.root","recreate");
 622
 623   for(Int_t iperiod=0;iperiod<nperiods;iperiod++){
 624     cout<<"Period "<<iperiod<<" : ";
 625     if(iperiod==0){ //LHC10bc
 626
 627       alephParameters[0] = 0.0283086/0.97;
 628       alephParameters[1] = 2.63394e+01;
 629       alephParameters[2] = 5.04114e-11;
 630       alephParameters[3] = 2.12543e+00;
 631       alephParameters[4] = 4.88663e+00;
 632       periodsname[0]="dataLHC10bc";
 633     }
 634     if(iperiod==1){ //LHC10de,low energy
 635       if(set==0){
 636         alephParameters[0] = 1.63246/50.;
 637         alephParameters[1] = 2.20028e+01;
 638         alephParameters[2] = TMath::Exp(-2.48879e+01);
 639         alephParameters[3] = 2.39804e+00;
 640         alephParameters[4] = 5.12090e+00;
 641         periodsname[1]="dataLHC10deold";
 642       }
 643       if(set==1){
 644         alephParameters[0] = 1.34490e+00/50.;
 645         alephParameters[1] =  2.69455e+01;
 646         alephParameters[2] =  TMath::Exp(-2.97552e+01);
 647         alephParameters[3] = 2.35339e+00;
 648         alephParameters[4] = 5.98079e+00;
 649         periodsname[2]="dataLHC10denew";
 650       }
 651     }
 652
 653     if(iperiod==2){//LHC10h
 654       if(set==0){//pass1
 655         alephParameters[0]=1.25202/50.;
 656         alephParameters[1]=2.74992e+01;
 657         alephParameters[2]=TMath::Exp(-3.31517e+01);
 658         alephParameters[3]=2.46246;
 659         alephParameters[4]=6.78938;
 660         periodsname[3]="dataLHC10hpass1";
 661       }
 662       if (set==1){//pass2 (AOD044)
 663         alephParameters[0]=1.25202/50.;
 664         alephParameters[1]=2.74992e+01;
 665         alephParameters[2]=TMath::Exp(-3.31517e+01);
 666         alephParameters[3]=2.46246;
 667         alephParameters[4]=6.78938;
 668         periodsname[4]="dataLHC10hpass2";
 669       }
 670     }
 671     if(iperiod==3){ //MC
 672       if(set==0){
 673         alephParameters[0] = 2.15898e+00/50.;
 674         alephParameters[1] = 1.75295e+01;
 675         alephParameters[2] = 3.40030e-09;
 676         alephParameters[3] = 1.96178e+00;
 677         alephParameters[4] = 3.91720e+00;
 678         periodsname[5]="MCold";
 679       }
 680       if(set==1){ //new
 681         alephParameters[0] = 1.44405/50;
 682         alephParameters[1] = 2.35409e+01;
 683         alephParameters[2] = TMath::Exp(-2.90330e+01);
 684         alephParameters[3] = 2.10681;
 685         alephParameters[4] = 4.62254;
 686         periodsname[6]="MCnew";
 687       }
 688
 689       if(set==2){ //new BB from Francesco
 690         alephParameters[0] = 0.029021;
 691         alephParameters[1] = 25.4181;
 692         alephParameters[2] = 4.66596e-08;
 693         alephParameters[3] = 1.90008;
 694         alephParameters[4] = 4.63783;
 695         periodsname[7]="MCBBFrancesco";
 696       }
 697
 698       if(set==3){ //low energy 2011
 699         alephParameters[0] = 0.0207667;
 700         alephParameters[1] = 29.9936;
 701         alephParameters[2] = 3.87866e-11;
 702         alephParameters[3] = 2.17291;
 703         alephParameters[4] = 7.1623;
 704         //periodsname[8]="MClowen2011";
 705       }
 706
 707
 708     }
 709     //cout<<periodsname[iperiod]<<endl;
 710     tpcResp->SetBetheBlochParameters(alephParameters[0],alephParameters[1],alephParameters[2],alephParameters[3],alephParameters[4]);
 711     cout<<"here"<<endl;
 712     for(Int_t ipart=0;ipart<npart;ipart++){
 713
 714       const Int_t n=1000;
 715       Double_t p[n],bethe[n];
 716
 717       for(Int_t k=0;k<n;k++){ //loop on the momentum steps
 718         p[k]=0.0001+k*4./n; //limits 0.-4. GeV/c
 719         //cout<<p[k]<<"\t";
 720         //bethe[k]=-tpcResp->Bethe(p[k]/masses[ipart]);
 721         AliPID::EParticleType ptype=AliPID::kKaon;
 722         if(ipart==1) ptype=AliPID::kPion;
 723         if(ipart==2) ptype=AliPID::kElectron;
 724         if(ipart==3) ptype=AliPID::kProton;
 725         bethe[k]=tpcResp->GetExpectedSignal(p[k],ptype);
 726       }
 727       //cout<<endl;
 728       TGraph *gr=new TGraph(n,p,bethe);
 729       gr->SetName(Form("%sP%d",partnames[ipart].Data(),iperiod));
 730       gr->SetTitle(Form("%sP%d;p (GeV/c);",partnames[ipart].Data(),iperiod));
 731       gr->SetLineColor(ipart+1);
 732       gr->SetMarkerColor(ipart+1);
 733       gr->GetYaxis()->SetRangeUser(35,100);
 734       cBethe->cd();
 735       if(iperiod==0 && ipart==0)gr->DrawClone("AL");
 736       else gr->DrawClone("L");
 737
 738       fout->cd();
 739       gr->Write();
 740     }
 741
 742   }
 743   TParameter<int> sett;
 744   sett.SetVal(set);
 745   fout->cd();
 746   sett.Write();
 747
 748   fout->Close();
 749 }
 750
 751 void DrawOutputCentrality(TString partname="D0",TString textleg="",TString path="./", Bool_t superimpose=kFALSE,TString suffixdir="",TString filename=/*"AnalysisResults.root"*/"PWG3histograms.root"){
 752   gStyle->SetCanvasColor(0);
 753   gStyle->SetTitleFillColor(0);
 754   gStyle->SetStatColor(0);
 755   gStyle->SetPalette(1);
 756
 757   TString listname="outputCentrCheck",partname2="",path2="",suffixdir2="",filename2="PWG3histograms.root";
 758
 759   // if(superimpose){
 760   //   cout<<"Enter the names:\n>";
 761   //   cin>>name1;
 762   //   cout<<">";
 763   //   cin>>name2;
 764   // }
 765   // TString listname="outputTrack",name1="",name2="";
 766   TString tmp="y";
 767
 768   if(superimpose){
 769     cout<<"##Second file\n";
 770     cout<<"Enter the name:\n";
 771     cout<<">";
 772     cin>>partname2;
 773     cout<<"Are they in the same output file? (y/n)"<<endl;
 774     cin>>tmp;
 775     if(tmp=="n"){
 776       cout<<"Path: \n";
 777       cout<<">";
 778       cin>>path2;
 779       cout<<"Dir name:\n";
 780       cout<<">";
 781       cin>>suffixdir2;
 782       cout<<"Filename: "<<endl;
 783       cout<<">";
 784       cin>>filename2;
 785     }
 786
 787   }
 788   // Int_t nhist=1;
 789   // TString *name=0x0;
 790   // if(superimpose){
 791   //   cout<<"Number of histogram to superimpose: ";
 792   //   cin>>nhist;
 793   //   name=new TString[nhist];
 794   //   for (Int_t j=0;j<nhist;j++){
 795   //     cout<<">";
 796   //     cin>>name[j];
 797   //   }
 798   // }
 799
 800   TList* list;
 801   TH1F * hstat;
 802
 803   TString dirname="PWG3_D2H_QA",dirname2=dirname;
 804   dirname+=suffixdir;
 805   dirname2+=suffixdir2;
 806   Bool_t isRead=ReadFile(list,hstat,listname,partname.Data(),path,filename,dirname);
 807   if(!isRead) return;
 808   if(!list || !hstat){
 809     cout<<":-( null pointers..."<<endl;
 810     return;
 811   }
 812
 813   TPaveText *pvtxt=new TPaveText(0.6,0.6,0.9,0.9,"NDC");
 814   pvtxt->SetBorderSize(0);
 815   pvtxt->SetFillStyle(0);
 816   pvtxt->AddText(partname);
 817
 818   TList* llist;
 819   TH1F* hhstat;
 820   if(superimpose){
 821     isRead=ReadFile(llist,hhstat,listname,partname2.Data(),path2,filename2,dirname2);
 822     if(!isRead) return;
 823     if(!llist || !hhstat){
 824       cout<<":-( null pointers..."<<endl;
 825       return;
 826     }
 827     TText *redtext=pvtxt->AddText(partname2);
 828     redtext->SetTextColor(kRed);
 829
 830   }
 831
 832
 833   TCanvas* cst=new TCanvas("cst","Stat");
 834   cst->SetGridy();
 835   cst->cd();
 836   Int_t nevents=hstat->GetBinContent(1);
 837   hstat->Draw("htext0");
 838   cst->SaveAs(Form("%s%s.png",hstat->GetName(),textleg.Data()));
 839   cst->SaveAs(Form("%s%s.eps",hstat->GetName(),textleg.Data()));
 840   Int_t nevents080=1,nnevents080=1;
 841
 842   //TCanvas *spare=new TCanvas("sparecv","Spare");
 843
 844   for(Int_t i=0;i<list->GetEntries();i++){
 845
 846     TClass* objtype=list->At(i)->IsA();
 847     TString tpname=objtype->GetName();
 848
 849     if(tpname=="TH1F"){
 850
 851       TH1F* h=(TH1F*)list->At(i);
 852       TH1F* hh=0x0;
 853       if(superimpose){
 854         hh=(TH1F*)llist->At(i);
 855       }
 856       if(!h || (superimpose && !hh)){
 857         cout<<"Histogram "<<i<<" not found"<<endl;
 858         continue;
 859       }
 860       if(superimpose){
 861         hhstat->SetLineColor(kRed);
 862         hh->SetLineColor(kRed);
 863       }
 864
 865       TCanvas* c=new TCanvas(Form("c%s",h->GetName()),h->GetName());
 866       TPaveText *pvtxt2=new TPaveText(0.6,0.6,0.9,0.9,"NDC");
 867       pvtxt2->SetBorderSize(0);
 868       pvtxt2->SetFillStyle(0);
 869
 870       c->cd();
 871       c->SetGrid();
 872       c->SetLogy();
 873       Int_t entries=h->Integral();
 874       pvtxt2->AddText(Form("%.1f %s of the events",(Double_t)entries/(Double_t)nevents*100,"%"));
 875       h->Draw();
 876       if(superimpose) {
 877         hh->Draw("sames");
 878         pvtxt->Draw();
 879       }
 880       pvtxt2->Draw();
 881       c->SaveAs(Form("%s%s.pdf",c->GetName(),textleg.Data()));
 882       c->SaveAs(Form("%s%s.eps",c->GetName(),textleg.Data()));
 883     }
 884     if(tpname=="TH2F"){
 885       TH2F* hhh=(TH2F*)list->At(i);
 886       if(!hhh){
 887         cout<<"Histogram "<<i<<" not found"<<endl;
 888         continue;
 889       }
 890       TCanvas* c=new TCanvas(Form("c%s",hhh->GetName()),hhh->GetName());
 891       TPaveText *pvtxt3=new TPaveText(0.6,0.6,0.9,0.9,"NDC");
 892       pvtxt3->SetBorderSize(0);
 893       pvtxt3->SetFillStyle(0);
 894
 895       c->cd();
 896       c->SetGrid();
 897       Int_t entries=hhh->Integral();
 898       pvtxt3->AddText(Form("%.1f %s of the events",(Double_t)entries/(Double_t)nevents*100,"%"));
 899       hhh->Draw("colz");
 900       c->SetLogz();
 901       pvtxt3->Draw();
 902       c->SaveAs(Form("%s%s.pdf",c->GetName(),textleg.Data()));
 903       c->SaveAs(Form("%s%s.eps",c->GetName(),textleg.Data()));
 904     }
 905   }
 906
 907
 908   listname="countersCentrality";
 909
 910   isRead=ReadFile(list,hstat,listname,partname.Data(),path,filename,dirname);
 911   if(!isRead) return;
 912   if(!list || !hstat){
 913     cout<<":-( null pointers..."<<endl;
 914     return;
 915   }
 916
 917
 918   if(superimpose){
 919     isRead=ReadFile(llist,hhstat,listname,partname2.Data(),path2,filename2,dirname2);
 920     if(!isRead) return;
 921     if(!llist || !hhstat){
 922       cout<<":-( null pointers..."<<endl;
 923       return;
 924     }
 925     TText *redtext=pvtxt->AddText(partname2);
 926     redtext->SetTextColor(kRed);
 927
 928   }
 929
 930   TH1F* hallcntr=0x0;
 931   TH1F* hhallcntr=0x0;
 932   cout<<"normalizing to 0-80% as a check"<<endl;
 933   Int_t ncentr=10;//check this
 934   TH1F* h020=0x0;
 935   TH1F* h2080=0x0;
 936   TH1F* hh020=0x0;
 937   TH1F* hh2080=0x0;
 938
 939   TCanvas *cvnocnt=new TCanvas("cvnocnt","No Centrality estimation",800,400);
 940   cvnocnt->Divide(2,1);
 941   TCanvas *ccent=0x0;
 942
 943   for(Int_t i=0;i<list->GetEntries();i++){
 944     AliCounterCollection* coll=(AliCounterCollection*)list->At(i);
 945     AliCounterCollection* colle=0x0;
 946     if(superimpose) colle=(AliCounterCollection*)llist->At(i);
 947     coll->SortRubric("run");//sort by run number
 948
 949     h020=0x0;
 950     h2080=0x0;
 951     hh020=0x0;
 952     hh2080=0x0;
 953
 954     hallcntr=0x0;
 955     hhallcntr=0x0;
 956
 957     TH1F* hbad=(TH1F*)coll->Get("run",Form("centralityclass:-990_-980"));
 958     cvnocnt->cd(i+1);
 959     if(hbad) hbad->Draw();
 960
 961     ccent=new TCanvas(Form("ccent%s",coll->GetName()),Form("Centrality vs Run (%s)",coll->GetName()),1400,800);
 962     ccent->SetTicky();
 963     ccent->Divide(4,2);
 964
 965     TH1F* hh=0x0;
 966
 967     for(Int_t ic=0;ic<8/*ncentr*/;ic++){ //normalizing to 0-80% as a check
 968
 969       TH1F* h=(TH1F*)coll->Get("run",Form("centralityclass:%d_%d",ic*10,ic*10+10));
 970       h->SetName(Form("h%d%d",i,ic));
 971       if(!hallcntr) {
 972         hallcntr=(TH1F*)h->Clone("hallcntr");
 973         hallcntr->Sumw2();
 974       } else {
 975         hallcntr->Add(h);
 976       }
 977
 978       nevents080+=h->Integral();
 979
 980       if(superimpose){
 981         hh=(TH1F*)colle->Get("run",Form("centralityclass:%d_%d",ic*10,ic*10+10));
 982         hh->SetName(Form("hh%d%d",i,ic));
 983         if(!hhallcntr) {
 984           hhallcntr=(TH1F*)hh->Clone("hhallcntr");
 985           hhallcntr->Sumw2();
 986         }else hhallcntr->Add(hh);
 987
 988         nnevents080+=hh->Integral();
 989
 990       }
 991     }
 992
 993     for(Int_t ic=0;ic<ncentr;ic++){
 994
 995       TH1F* h=(TH1F*)coll->Get("run",Form("centralityclass:%d_%d",ic*10,ic*10+10));
 996       h->SetName(Form("h%d%d",i,ic));
 997       h->Sumw2();
 998
 999       if(ic>=0 && ic<=1){ //0-20
1000         if(!h020) {
1001           h020=(TH1F*)h->Clone(Form("h020%s",coll->GetName()));
1002           h020->SetTitle(Form("Centrality 0-20 %s",coll->GetName()));
1003           if(superimpose){
1004             hh020=(TH1F*)hh->Clone(Form("hh020%s",coll->GetName()));
1005             hh020->SetTitle(Form("Centrality 0-20 %s",coll->GetName()));
1006           }
1007         }
1008         else {
1009           h020->Add(h);
1010           if(superimpose)hh020->Add(hh);
1011         }
1012       }
1013       if(ic>=2 && ic<=7){ //20-80
1014         if(!h2080) {
1015           h2080=(TH1F*)h->Clone(Form("h2080%s",coll->GetName()));
1016           h2080->SetTitle(Form("Centrality 20-80 %s",coll->GetName()));
1017           if(superimpose){
1018             hh2080=(TH1F*)hh->Clone(Form("hh2080%s",coll->GetName()));
1019             hh2080->SetTitle(Form("Centrality 20-80 %s",coll->GetName()));
1020           }
1021         }
1022         else {
1023           h2080->Add(h);
1024           if(superimpose)hh2080->Add(hh);
1025         }
1026
1027       }
1028
1029       h->Divide(hallcntr);
1030
1031       if(ic<8){
1032         ccent->cd(ic+1);
1033         h->GetYaxis()->SetLabelSize(0.05);
1034         h->GetYaxis()->SetTitleOffset(1.5);
1035         h->SetMinimum(0);
1036         //h->GetYaxis()->SetRangeUser(0.,0.15);
1037         h->DrawClone();
1038       }
1039       /*
1040         if(ic==0&&i==0){
1041         spare->cd();
1042         h->Draw();
1043         }
1044       */
1045       // ccent->cd(1);
1046       // h->SetLineColor(ic+1);
1047       // if(ic==0)h->DrawClone();
1048       // else h->DrawClone("sames");
1049     }
1050     h020->Divide(hallcntr);
1051     if(superimpose){
1052       hh020->Divide(hhallcntr);
1053       hh020->SetLineColor(2);
1054       hh020->SetMarkerColor(2);
1055     }
1056
1057     /*//draw 0-20 and 20-80 in the multi pad canvas (increase divisions before uncommenting)
1058     ccent->cd(ncentr+1);
1059     h020->DrawClone();
1060     if(superimpose){
1061       hh020->DrawClone("sames");
1062     }
1063     */
1064     TCanvas* cv020=new TCanvas(Form("cv020-%d",i),"0-20% vs run number",1400,600);
1065     cv020->cd();
1066     cv020->SetTicky();
1067     h020->GetYaxis()->SetRangeUser(0.,1.);
1068     h020->DrawClone();
1069     if(superimpose)hh020->DrawClone("sames");
1070     cv020->SaveAs(Form("cv020-%d.pdf",i));
1071     cv020->SaveAs(Form("cv020-%d.eps",i));
1072
1073     h2080->Divide(hallcntr);
1074     if(superimpose) {
1075       hh2080->Divide(hhallcntr);
1076       hh2080->SetLineColor(2);
1077       hh2080->SetMarkerColor(2);
1078     }
1079
1080     /*
1081     ccent->cd(ncentr+2);
1082     h2080->DrawClone();
1083
1084     if(superimpose){
1085       hh2080->DrawClone("sames");
1086     }
1087     */
1088     TCanvas* cv2080=new TCanvas(Form("cv2080-%d",i),"20-80% vs run number",1400,600);
1089     cv2080->cd();
1090     cv2080->SetTicky();
1091     h2080->GetYaxis()->SetRangeUser(0.,1.);
1092     h2080->DrawClone();
1093     if(superimpose)hh2080->DrawClone("sames");
1094     cv2080->SaveAs(Form("cv2080-%d.pdf",i));
1095     cv2080->SaveAs(Form("cv2080-%d.eps",i));
1096
1097     ccent->SaveAs(Form("%s%s.pdf",ccent->GetName(),textleg.Data()));
1098     ccent->SaveAs(Form("%s%s.eps",ccent->GetName(),textleg.Data()));
1099   }
1100
1101 }
1102
1103 void DrawProjections(TString partname="D0",TString h2dname="hMultvsPercentile",Int_t groupnbins=5,Float_t fitmin=15,Float_t fitmax=50,TString direction="X",TString path="./",TString suffixdir="", TString filename="AnalysisResults.root", TString fitfunc="pol0"/*option "nofit" does not fit*/){
1104   gStyle->SetCanvasColor(0);
1105   gStyle->SetTitleFillColor(0);
1106   gStyle->SetStatColor(0);
1107   gStyle->SetPalette(1);
1108
1109   TString listname="outputCentrCheck";
1110   TString dirname="PWG3_D2H_QA";
1111   dirname+=suffixdir;
1112
1113   TList* list;
1114   TH1F * hstat;
1115
1116   Bool_t isRead=ReadFile(list,hstat,listname,partname,path,filename,dirname);
1117   if(!isRead) return;
1118   if(!list || !hstat){
1119     cout<<":-( null pointers..."<<endl;
1120     return;
1121   }
1122   Double_t nevents=hstat->GetBinContent(5); //ev good vertex
1123
1124   TH2F* h2=(TH2F*)list->FindObject(h2dname);
1125   if(!h2){
1126     cout<<h2dname.Data()<<" not found"<<endl;
1127     list->ls();
1128     return;
1129   }
1130   TCanvas* cv2d=new TCanvas("cv2d",h2->GetName());
1131   cv2d->cd();
1132   cv2d->SetLogz();
1133   cv2d->SetGrid();
1134   h2->Draw("colz");
1135   TPaveText *pvst=new TPaveText(0.6,0.2,0.9,0.7,"NDC");
1136   pvst->SetBorderSize(0);
1137   pvst->SetFillStyle(0);
1138   pvst->AddText("Bin -> Cont/nEvVtx");
1139
1140   //nsteps=group bins in the Y(X) direction if projecting on the X(Y) direction
1141   Int_t nsteps=0;
1142
1143   if(direction=="X") nsteps=h2->GetNbinsY()/groupnbins;
1144   if(direction=="Y") nsteps=h2->GetNbinsX()/groupnbins;
1145   cout<<"Grouping bins by " <<groupnbins<<" I obtaine "<<nsteps<<" projections"<<endl;
1146
1147   TCanvas *cvpj=new TCanvas(Form("cvpj%s%s",direction.Data(),h2dname.Data()),Form("cvpj%s",direction.Data()),1200,800);
1148   cvpj->Divide((Int_t)(nsteps/3)+1,3);
1149   TFile* fout=new TFile(Form("proj%s%s.root",direction.Data(),h2dname.Data()), "recreate");
1150   //Float_t maxx[nsteps];
1151   //Float_t maxx[12]={9000,9000,6000,4000,2000,1400,800,500,200,100,40,25};
1152   Double_t integralpernev[nsteps];
1153
1154   Double_t minx=0,maxx=0;
1155   if(direction=="X"){
1156     minx=h2->GetYaxis()->GetXmin();
1157     maxx=h2->GetYaxis()->GetXmax();
1158   }
1159   if(direction=="Y"){
1160     minx=h2->GetXaxis()->GetXmin();
1161     maxx=h2->GetXaxis()->GetXmax();
1162   }
1163   printf("Plotting from %.1f to %.1f\n",minx,maxx);
1164   TCanvas *cintegral=new TCanvas("cintegral","Integral of each projection");
1165   TH1F* hint=new TH1F("hint","Integral of each projection;Centrality (%);Entries",nsteps,minx,maxx);
1166   Double_t minint=999999999,maxint=0;
1167
1168   for(Int_t i=0;i<nsteps;i++){
1169     TH1F* h=0x0;
1170     // if(direction=="X")h=(TH1F*)h2->ProjectionX(Form("px%d",i),i+kbins,i+2*kbins);
1171     // if(direction=="Y")h=(TH1F*)h2->ProjectionY(Form("py%d",i),i+kbins,i+2*kbins);
1172     if(direction=="X")h=(TH1F*)h2->ProjectionX(Form("px%d",i),groupnbins*i+1,groupnbins*(i+1));
1173     if(direction=="Y")h=(TH1F*)h2->ProjectionY(Form("py%d",i),groupnbins*i+1,groupnbins*(i+1));
1174     Double_t projint=h->Integral();
1175     cout<<"Integral of projection "<<i<<" = "<<projint<<endl;
1176     hint->SetBinContent(i+1,projint);
1177     hint->SetBinError(i+1,TMath::Sqrt(projint));
1178
1179     if(projint<1e-7) continue;
1180     if(minint>projint) minint=projint;
1181     if(projint>maxint) maxint=projint;
1182     integralpernev[i]=h->Integral()/nevents;
1183
1184     TPaveText *pvtxt=new TPaveText(0.6,0.6,0.9,0.9,"NDC");
1185     pvtxt->SetBorderSize(0);
1186     pvtxt->SetFillStyle(0);
1187     pvtxt->AddText(Form("%.0f - %.0f",h2->GetYaxis()->GetBinLowEdge(groupnbins*i+1),h2->GetYaxis()->GetBinLowEdge(groupnbins*(i+1))));
1188     pvst->AddText(Form("%.0f - %.0f -> %.2f",h2->GetYaxis()->GetBinLowEdge(groupnbins*i+1),h2->GetYaxis()->GetBinLowEdge((groupnbins*(i+1))),integralpernev[i]));
1189
1190     cvpj->cd(i+1);
1191     //h->GetXaxis()->SetRangeUser(0,maxx[i]);
1192     h->Draw();
1193     pvtxt->Draw();
1194     fout->cd();
1195     h->Write();
1196   }
1197   cvpj->SaveAs(Form("cvpj%s%s.pdf",direction.Data(),h2dname.Data()));
1198   cvpj->SaveAs(Form("cvpj%s%s.eps",direction.Data(),h2dname.Data()));
1199
1200   cv2d->cd();
1201   pvst->Draw();
1202   cv2d->SaveAs(Form("%s.pdf",h2->GetName()));
1203   cv2d->SaveAs(Form("%s.eps",h2->GetName()));
1204
1205   cintegral->cd();
1206   hint->SetMarkerStyle(20);
1207   hint->Draw("PE");
1208   if(!fitfunc.Contains("nofit")){
1209     hint->Fit(fitfunc.Data(),"RL","PE",fitmin,fitmax);
1210     TF1* fpolfit=hint->GetFunction(fitfunc.Data());
1211     TPaveText *txtvar=new TPaveText(0.3,0.1,0.9,0.4,"NDC");
1212     txtvar->SetBorderSize(0);
1213     txtvar->SetFillStyle(0);
1214     //txtvar->AddText(Form("Full spread %.0f- %.0f",maxint,minint));
1215     txtvar->AddText(Form("Fit in %.1f-%.1f; ",fitmin,fitmax));
1216     for(Int_t ipar=0;ipar<fpolfit->GetNpar();ipar++){
1217       txtvar->AddText(Form("par%d = %.0f, ",ipar, fpolfit->GetParameter(ipar)));
1218     }
1219     txtvar->AddText(Form("#tilde{#chi}^{2} = %.2f",fpolfit->GetChisquare()/fpolfit->GetNDF()));
1220     txtvar->AddText(Form("bin width = %.1f %s",hint->GetBinWidth(3),"%"));
1221     txtvar->Draw();
1222   }
1223   fout->cd();
1224   hint->Write();
1225   cintegral->SaveAs(Form("%s.pdf",hint->GetName()));
1226   cintegral->SaveAs(Form("%s.eps",hint->GetName()));
1227 }
1228
1229 void DrawEventSelection(TString partname="D0", TString path="./",TString suffixdir="",TString filename="AnalysisResults.root"){
1230   gStyle->SetCanvasColor(0);
1231   gStyle->SetTitleFillColor(0);
1232   gStyle->SetStatColor(0);
1233   gStyle->SetPalette(1);
1234   gStyle->SetOptStat(0);
1235
1236   TString listname="outputEvSel";
1237   TString dirname="PWG3_D2H_QA";
1238   dirname+=suffixdir;
1239
1240   TList* list;
1241   TH1F * hstat;
1242
1243   Bool_t isRead=ReadFile(list,hstat,listname,partname,path,filename,dirname);
1244   if(!isRead) return;
1245   if(!list || !hstat){
1246     cout<<":-( null pointers..."<<endl;
1247     return;
1248   }
1249   //Double_t neventsgv=hstat->Integral(5,5); //ev good vertex
1250
1251   for(Int_t i=0;i<list->GetEntries();i++){
1252
1253     TClass* objtype=list->At(i)->IsA();
1254     TString tpname=objtype->GetName();
1255
1256     if(tpname=="TH1F"){
1257       TH1F* htmp=(TH1F*)list->At(i);
1258       TCanvas* c=new TCanvas(Form("c%s",htmp->GetName()),Form("c%s",htmp->GetName()));
1259       c->cd();
1260       htmp->Draw();
1261       c->SaveAs(Form("%s.pdf",htmp->GetName()));
1262       c->SaveAs(Form("%s.eps",htmp->GetName()));
1263     }
1264
1265     if(tpname=="TH2F"){
1266       TH2F* htmp=(TH2F*)list->At(i);
1267       TCanvas* c=new TCanvas(Form("c%s",htmp->GetName()),Form("c%s",htmp->GetName()),1200,800);
1268       c->cd();
1269       htmp->SetMarkerSize(1.3);
1270       htmp->Draw("colzhtext45");
1271       c->SaveAs(Form("%s.pdf",htmp->GetName()));
1272       c->SaveAs(Form("%s.eps",htmp->GetName()));
1273     }
1274   }
1275
1276   AliCounterCollection* coll=(AliCounterCollection*)list->FindObject("trigCounter");
1277   if(!coll) {
1278     cout<<"Trigger counter not found"<<endl;
1279     return;
1280   }
1281
1282   coll->SortRubric("run");//sort by run number
1283
1284   TString collname=coll->GetName();
1285
1286   TString keywords=coll->GetKeyWords("triggertype");
1287
1288   Int_t nkeys=keywords.CountChar(',')+1;
1289
1290   TString *words = new TString[nkeys];
1291   for(Int_t k=0;k<nkeys;k++) words[k]="";
1292   printf("Keywords: ");
1293   Int_t count=0;
1294   for(Int_t l=0;l<keywords.Length();l++){
1295     if(keywords[l] != ',') words[count]+=keywords[l];
1296     else {
1297       printf("%s ",words[count].Data());
1298       count++;
1299     }
1300   }
1301   cout<<endl;
1302
1303   TH1D** htrig=new TH1D*[nkeys]; //each trigger type in one histogram of counts vs run
1304   TH1D** htrignorm=new TH1D*[nkeys]; //normalized to the counts in kAny
1305   TCanvas* ctrigfraction=new TCanvas("cvtrigfrac","Fraction of given trigger type vs run",1400,800);
1306   TLegend* legtr=new TLegend(0.15,0.5,0.35,0.8);
1307   legtr->SetBorderSize(0);
1308   legtr->SetFillStyle(0);
1309   for(Int_t k=0;k<nkeys;k++){
1310     htrig[k]=coll->Get("run",Form("triggertype:%s",words[k].Data()));
1311     htrig[k]->SetName(Form("h%s",words[k].Data()));
1312     htrig[k]->SetTitle("Trigger type;RUN; counts");
1313     htrig[k]->SetMarkerColor(k+1);
1314     htrig[k]->SetMarkerStyle(k+20);
1315     htrig[k]->Sumw2();
1316     legtr->AddEntry(htrig[k],Form("%s",words[k].Data()),"P");
1317     //drawings
1318     //1) counts of a given trigger over counts in kAny
1319     htrignorm[k]=(TH1D*)htrig[k]->Clone(Form("h%snormAny",words[k].Data()));
1320     htrignorm[k]->SetTitle("Trigger type over ANY trigger;RUN; counts/countsANY");
1321     htrignorm[k]->Divide(htrig[k],htrig[0],1.,1.,"B");
1322     htrignorm[k]->GetXaxis()->SetRangeUser(0,1.1);
1323
1324     ctrigfraction->cd();
1325     if(k>0)htrignorm[k]->Draw("PEsames");
1326     else htrignorm[k]->Draw("PE");
1327   }
1328
1329   ctrigfraction->cd();
1330   legtr->Draw();
1331   ctrigfraction->SaveAs("TrgFractionOverANY.pdf");
1332   ctrigfraction->SaveAs("TrgFractionOverANY.eps");
1333
1334   delete [] words;
1335
1336   //draw number of events per trigger
1337   TH2F* hTrigMul=(TH2F*)list->FindObject("hTrigMul");
1338
1339   TH1F *hnEvPerTrig=(TH1F*)hTrigMul->ProjectionX("hnEvPerTrig");
1340   hnEvPerTrig->SetTitle("Number of events per trigger");
1341   TCanvas* cnev=new TCanvas("cnev","Number of events",1400,800);
1342   cnev->cd();
1343   hnEvPerTrig->Draw("htext0");
1344   cnev->SaveAs(Form("%s.pdf",cnev->GetName()));
1345   cnev->SaveAs(Form("%s.eps",cnev->GetName()));
1346
1347   //draw number of events selected per trigger
1348   TH2F* hTrigMulSel=(TH2F*)list->FindObject("hTrigMulSel");
1349
1350   TH1F *hnEvPerTrigSel=(TH1F*)hTrigMulSel->ProjectionX("hnEvPerTrigSel");
1351   hnEvPerTrigSel->SetTitle("Number of events selected per trigger");
1352   TCanvas* cnevs=new TCanvas("cnevs","Number of events selected",1400,800);
1353   cnevs->cd();
1354   hnEvPerTrigSel->Draw("htext0");
1355   cnevs->SaveAs(Form("%s%s.pdf",cnevs->GetName(),suffixdir.Data()));
1356   cnevs->SaveAs(Form("%s%s.eps",cnevs->GetName(),suffixdir.Data()));
1357
1358   TH1F* hnEvSelPerc=(TH1F*)hnEvPerTrigSel->Clone("hnEvSelFrac");
1359   hnEvSelPerc->SetTitle("Fraction of event selected per trigger");
1360   hnEvSelPerc->Divide(hnEvPerTrig);
1361   TCanvas* cnevsp=new TCanvas("cnevsp","Fraction of events selected",1400,800);
1362   cnevsp->cd();
1363   hnEvSelPerc->Draw("htext0");
1364   cnevsp->SaveAs(Form("%s%s.pdf",cnevsp->GetName(),suffixdir.Data()));
1365   cnevsp->SaveAs(Form("%s%s.eps",cnevsp->GetName(),suffixdir.Data()));
1366 }
1367
1368 void WriteTextFileForCompilation(TString inputtextfile, Int_t analysistype){
1369    //This method writes a text file with the name given in input containing the path of the QA output file, the directories and the lists of histograms to be read to perform a comparison of QA outputs
1370    // analysistype = 1 is the OutputTrack
1371    // analysistype = 2 is the OutputEventSelection
1372    // other outputs can be implemented
1373
1374    // The paths and trigger type (or more in general the suffix of the directories to be read have to be set in this method)
1375    // The legend is set either equal to the trigger type if they are different, or it should be set manually
1376    // Run this methos before CompilationOfTriggeredEvents unless the text file has been already produced
1377
1378   Printf("WARNING! Did you customize the parameters in the macro?");
1379   const Int_t n=8; //customize this
1380   TString prepath="/data/Work/HFQA/LHC12QATrainOutputs/";//"./";
1381   TString pathname[n]={"LHC12a/Train165/","LHC12b/Train166/","LHC12d/Train168/","LHC12e/Train170/","LHC12f/Train171/","LHC12g/Train172/"};
1382
1383   TString dirname[n];
1384   TString listname[n];
1385   ofstream myfile;
1386   myfile.open(inputtextfile);
1387   myfile<<n<<endl;
1388
1389   TString filename="AnalysisResults.root";
1390   TString basedirname="PWG3_D2H_QA";
1391   TString trgnames[n];//={"INT7","EMC7","EMCJET7","EMCGAMMA7","INT8","EMC8","EMCJET8","EMCGAMMA8"}; //customize this
1392   TString baselistname="", partname="D00100";
1393   TString legend[n]={"LHC12a165","LHC12b166","LHC12d168","LHC12e170","LHC12f171","LHC12g172"}; //Customize this if want it to be different from trgnames
1394   if(analysistype==1) baselistname="outputTrack";
1395   if(analysistype==2) baselistname="outputEvSel";
1396   //... others could be added
1397
1398   //All in the same directory on your computer
1399   for(Int_t i=0;i<n;i++){
1400     //set
1401     pathname[i]=prepath+pathname[i];//"./";
1402     trgnames[i]="EMC7";
1403     dirname[i]=basedirname+trgnames[i];
1404     listname[i]=baselistname+partname+trgnames[i];
1405     if(legend[i]=="") legend[i]=trgnames[i];
1406
1407     //Printf("Trigger name is %s",trgnames[i].Data());
1408     //Printf("Legend is %s",legend[i].Data());
1409
1410     //write text file
1411     myfile<<endl;
1412     myfile<<endl;
1413     myfile<<pathname[i].Data()<<filename.Data()<<endl;
1414     myfile<<dirname[i].Data()<<endl;
1415     myfile<<listname[i].Data()<<endl;
1416     myfile<<legend[i].Data()<<endl;
1417
1418   }
1419   myfile.close();
1420
1421 }
1422
1423 Bool_t ReadFilesForCompilation(TString inputtextfile, TList**& lists, Int_t& numb, TString*& legend){
1424
1425    //Method to read the QA files indicated in the text file in input (to be written with WriteTextFileForCompilation() method)
1426   ifstream myfile;
1427   myfile.open(inputtextfile);
1428   if(!myfile.is_open()) Printf("No files found, did you make it correctly?");
1429   Int_t n;
1430   myfile >> n;
1431   lists=new TList*[n];
1432   legend= new TString[n];
1433   numb=n;
1434
1435   Int_t k=0;
1436
1437   while(myfile){
1438
1439     TString filename;
1440     myfile>>filename;
1441     //Printf("Filename = %s",filename.Data());
1442     TFile *fin=new TFile(filename);
1443     if(!fin->IsOpen()){
1444       Printf("File %s not found",filename.Data());
1445       return kFALSE;
1446     }
1447
1448     TString dirname;
1449     myfile>>dirname;
1450     //Printf("Dirname = %s",dirname.Data());
1451     TDirectoryFile* dir=(TDirectoryFile*)fin->Get(dirname);
1452     if(!dir){
1453       Printf("Directory %s not found in %s",dirname.Data(),filename.Data());
1454       fin->ls();
1455       return kFALSE;
1456     }
1457
1458     TString listname;
1459     myfile>>listname;
1460     //Printf("Listname = %s",listname.Data());
1461     lists[k]=(TList*)dir->Get(listname);
1462     if(!lists[k]){
1463       Printf("List %s not found in %s",listname.Data(),dirname.Data());
1464       dir->ls();
1465       return kFALSE;
1466     }
1467
1468     TString temp;
1469     myfile>>temp;
1470     legend[k]=temp;
1471     Printf("Legend = %s",legend[k].Data());
1472     if(!legend[k]){
1473       Printf("Legend %s not found",filename.Data());
1474       return kFALSE;
1475     }
1476
1477
1478     k++;
1479     if(k==n)return kTRUE;// Needed as the while loop does not always realise the end is reached.
1480   }
1481   return kTRUE;
1482 }
1483
1484 void CompilationOfTriggeredEvents(TString inputtextfile,Int_t analysistype){
1485    //This method draws the histograms superimposed
1486
1487    // analysistype = 1 is the OutputTrack
1488    // analysistype = 2 is the OutputEventSelection
1489    // other outputs can be implemented
1490
1491   TList **lists=0x0;
1492   Int_t n;
1493   TString* legend=0x0;
1494   Bool_t okread=ReadFilesForCompilation(inputtextfile,lists,n,legend);
1495   if(!okread){
1496     Printf("Did you write %s with  WriteTextFileForCompilation(%s)?",inputtextfile.Data(),inputtextfile.Data());
1497     return;
1498   }
1499
1500   if(analysistype==2){CompilationEvSelection(n,lists,legend);}
1501   if(analysistype==1){CompilationTrackSelection(n,lists,legend);}
1502 }
1503
1504 void CompilationEvSelection(Int_t n,TList** lists, TString* legend){
1505    // Specific method for EventSelection output (2)
1506
1507   TCanvas *cnevs=new TCanvas("cnevs","Number of events selected",1400,800);
1508   TCanvas *cnevsp=new TCanvas("cnevsp","Fraction of events selected",1400,800);
1509   TH1F** hnEvPerTrig=new TH1F*[n];
1510   TH1F** hnEvPerTrigSel=new TH1F*[n];
1511   TH1F** hnEvSelPerc=new TH1F*[n];
1512
1513   TLegend* leg=new TLegend(0.15,0.5,0.45,0.78);
1514   leg->SetFillStyle(0);
1515   leg->SetBorderSize(0);
1516
1517   Bool_t first=kTRUE;
1518   Int_t nentriesl=lists[0]->GetEntries();
1519   TCanvas** c=new TCanvas*[nentriesl];
1520   Int_t nth1f=0;
1521
1522   for(Int_t i=0;i<n;i++){
1523
1524     TH2F* hTrigMul=(TH2F*)lists[i]->FindObject("hTrigMul");
1525
1526     hnEvPerTrig[i]=(TH1F*)hTrigMul->ProjectionX(Form("hnEvPerTrig%d",i));
1527     hnEvPerTrig[i]->SetTitle("Number of events selected per trigger");
1528
1529
1530     TH2F* hTrigMulSel=(TH2F*)lists[i]->FindObject("hTrigMulSel");
1531
1532     hnEvPerTrigSel[i]=(TH1F*)hTrigMulSel->ProjectionX(Form("hnEvPerTrigSel%d",i));
1533     hnEvPerTrigSel[i]->SetTitle("Number of events selected per trigger");
1534     hnEvPerTrigSel[i]->SetLineColor(i+1);
1535     hnEvPerTrigSel[i]->SetLineWidth(2);
1536
1537     cnevs->cd();
1538     if(i==0) hnEvPerTrigSel[i]->Draw();
1539     else  hnEvPerTrigSel[i]->Draw("sames");
1540
1541     hnEvSelPerc[i]=(TH1F*)hnEvPerTrigSel[i]->Clone(Form("hnEvSelFrac%d",i));
1542     hnEvSelPerc[i]->SetTitle("Fraction of event selected per trigger");
1543     hnEvSelPerc[i]->SetLineColor(i+1);
1544     hnEvSelPerc[i]->SetLineWidth(2);
1545     hnEvSelPerc[i]->Divide(hnEvPerTrig[i]);
1546     cnevsp->cd();
1547     if(i==0) hnEvSelPerc[i]->Draw("htext0");
1548     else  hnEvSelPerc[i]->Draw("htext0sames");
1549     nth1f=0; //restart counting per each file
1550
1551     //fill legend
1552     leg->AddEntry(hnEvSelPerc[i],legend[i],"L");
1553
1554     for(Int_t j=0; j<nentriesl; j++){
1555        TClass* objtype=lists[i]->At(j)->IsA();
1556        TString tpname=objtype->GetName();
1557
1558        if (tpname=="TH1F"){
1559           TH1F* htmp=(TH1F*)lists[i]->At(j);
1560
1561           htmp->SetLineColor(1+i);
1562           if(htmp->Integral()>0) htmp->Scale(1./htmp->Integral());
1563
1564           if(first) {
1565              c[nth1f]=new TCanvas(Form("c%s",htmp->GetName()),Form("c%s",htmp->GetName()));
1566              c[nth1f]->cd();
1567              htmp->Draw();
1568
1569           }
1570           else {
1571              c[nth1f]->cd();
1572              htmp->Draw("sames");
1573           }
1574           nth1f++;
1575        }
1576     }
1577
1578     first=kFALSE;
1579
1580   }
1581
1582   for(Int_t j=0;j<nth1f;j++){
1583      c[j]->cd();
1584      leg->Draw();
1585      c[j]->SaveAs(Form("%s.pdf",c[j]->GetName()));
1586      c[j]->SaveAs(Form("%s.eps",c[j]->GetName()));
1587   }
1588
1589   cnevs->cd();
1590   leg->Draw();
1591   cnevs->SaveAs(Form("%s.eps",cnevs->GetName()));
1592   cnevs->SaveAs(Form("%s.pdf",cnevs->GetName()));
1593
1594   cnevsp->cd();
1595   leg->Draw();
1596   cnevsp->SaveAs(Form("%s.eps",cnevsp->GetName()));
1597   cnevsp->SaveAs(Form("%s.pdf",cnevsp->GetName()));
1598   TCanvas *test=new TCanvas();
1599   test->cd(); leg->Draw();
1600 }
1601
1602
1603 void CompilationTrackSelection(Int_t n,TList** lists,TString* legend)
1604 {
1605    // Specific method for Track output (1)
1606
1607   for(Int_t i=0;i<lists[0]->GetEntries();i++){
1608     TObjArray* plotseth=new TObjArray();
1609     TObjArray* plotsethr=new TObjArray();
1610
1611     Bool_t is2d=false;
1612     is2d = lists[0]->At(i)->InheritsFrom("TH2");
1613     if(is2d) continue;
1614
1615     for(Int_t j=0; j<n; j++){
1616       TH1F* temph=(TH1F*)lists[j]->At(i);
1617       TH1F * temphr=0x0;
1618       TH1F * ploto=0x0; //0th plot to divide by
1619
1620       // plot ratios and scale all plots to integral
1621       temphr= (TH1F*)temph->Clone(Form("%s_ratio",temph->GetName()));
1622       if(j==0)ploto = temphr;
1623       else ploto = (TH1F*)plotseth->At(0);
1624       temphr->Divide(ploto);
1625       if(temphr->Integral()>0){temphr-> Scale(1./temphr->Integral());}// take out for nonscaled ratios
1626       if(temph->Integral()>0){ temph-> Scale(1./temph->Integral());}
1627
1628       plotseth->AddLast(new TH1F(*temph));
1629       plotsethr->AddLast(new TH1F(*temphr));
1630     }
1631
1632
1633     TH1F *h = (TH1F*)plotseth->At(0);
1634     TCanvas* c=new TCanvas(Form("c%s",h->GetName()),h->GetName());
1635     c->cd();
1636     c->SetGrid();
1637     TString hname=h->GetName();
1638
1639     if(!hname.Contains("nCls")){
1640       c->SetLogy();
1641       if(hname.Contains("Layer")){
1642         for(Int_t ibin=1;ibin<=h->GetNbinsX();ibin++){
1643           h->GetXaxis()->SetBinLabel(ibin+1,Form("%d",ibin));
1644         }
1645         h->GetXaxis()->SetLabelSize(0.06);
1646         h->GetXaxis()->SetRangeUser(0,6); //comment to see ntracks!
1647       }
1648       //h->SetMinimum(1);
1649       h->Draw();
1650
1651       TLegend* leg=new TLegend(0.15,0.5,0.45,0.78);
1652       leg->SetFillStyle(0);
1653       leg->SetBorderSize(0);
1654       leg->AddEntry(h,legend[0],"L");
1655
1656       for(Int_t j=1; j<n; j++){
1657         TH1F * h2 = (TH1F*)plotseth->At(j);
1658         h2->SetLineColor(1+j);
1659         h2->Draw("sames");
1660         leg->AddEntry(h2,legend[j],"L");
1661       }
1662
1663       leg->Draw();
1664
1665       TCanvas* c2=new TCanvas(Form("c2%s",h->GetName()),h->GetName());
1666       c2->cd();
1667       c2->SetGrid();
1668       TH1F *hr = (TH1F*)plotsethr->At(1);
1669       hr->SetLineColor(2);
1670       hr->Draw();
1671       TLegend* leg2=new TLegend(0.15,0.5,0.45,0.78);
1672       leg2->SetFillStyle(0);
1673       leg2->SetBorderSize(0);
1674
1675       TString ratioleg;
1676       ratioleg+=legend[1];
1677       ratioleg+="/";
1678       ratioleg+=legend[0];
1679
1680       leg2->AddEntry(hr,ratioleg,"L");
1681
1682       for(Int_t j=2; j<n; j++){
1683         TH1F * hr2 = (TH1F*)plotsethr->At(j);
1684         hr2->SetLineColor(1+j);
1685         hr2->Draw("sames");
1686
1687         TString ratioleg2;
1688         ratioleg2+=legend[j];
1689         ratioleg2+="/";
1690         ratioleg2+=legend[0];
1691
1692         leg2->AddEntry(hr2,ratioleg2,"L");
1693       }
1694       leg2->Draw();
1695       c2->SaveAs(Form("%s%s%sRatio.pdf",c->GetName(),legend[0].Data(), legend[n-1].Data()));
1696       c2->SaveAs(Form("%s%s%sRatio.eps",c->GetName(),legend[0].Data(), legend[n-1].Data()));
1697     }
1698     else {
1699       h->Draw("htext0");
1700
1701       TLegend* leg=new TLegend(0.15,0.5,0.45,0.78);
1702       leg->SetFillStyle(0);
1703       leg->SetBorderSize(0);
1704       leg->AddEntry(h,legend[0],"L");
1705
1706       for(Int_t j=1; j<n; j++){
1707         TH1F * h2 = (TH1F*)plotseth->At(j);
1708         h2->SetLineColor(1+j);
1709         leg->AddEntry(h2,legend[j],"L");
1710         h2->Draw("htext0sames");
1711       }
1712
1713       leg->Draw();
1714     }
1715     c->cd();
1716     c->SaveAs(Form("%s%s%s.eps",c->GetName(),legend[0].Data(), legend[n-1].Data()));
1717     c->SaveAs(Form("%s%s%s.pdf",c->GetName(),legend[0].Data(), legend[n-1].Data()));
1718   }
1719
1720   Int_t check=0;
1721   Int_t maxfa=1;
1722   if(n<=10)maxfa=n;
1723   else cout<<"Warning: To many files for combinationplot, show only original"<<endl;
1724   TLegend* leg3=new TLegend(0.15,0.5,0.45,0.78);
1725   leg3->SetFillStyle(0);
1726   leg3->SetBorderSize(0);
1727   TCanvas* ctrsel=new TCanvas("ctrsel","Track Sel");
1728   ctrsel->SetLogy();
1729   for(Int_t i=0; i<maxfa; i++){
1730     TH1F* hd0fb4=(TH1F*)lists[i]->FindObject("hd0TracksFilterBit4");
1731     TH1F* hd0SPD1=(TH1F*)lists[i]->FindObject("hd0TracksSPDin");
1732     TH1F* hd0SPDany=(TH1F*)lists[i]->FindObject("hd0TracksSPDany");
1733     TH1F* hd0TPCITScuts=(TH1F*)lists[i]->FindObject("hd0TracksTPCITSSPDany");
1734     if(hd0fb4 && hd0SPD1 && hd0SPDany && hd0TPCITScuts){
1735       if(i==0){check=1;}
1736       else{if(check==0)return;}
1737       ctrsel->cd();
1738       hd0SPD1->SetLineColor(kCyan+3);
1739       hd0SPDany->SetLineColor(4);
1740       hd0TPCITScuts->SetLineColor(kGreen+1);
1741       hd0fb4->SetLineColor(2);
1742       if(i==0){
1743         hd0SPD1->Draw();
1744         ctrsel->Update();
1745         TPaveStats *st1=(TPaveStats*)hd0SPD1->GetListOfFunctions()->FindObject("stats");
1746         st1->SetTextColor(kCyan+3);
1747         st1->SetY1NDC(0.71);
1748         st1->SetY2NDC(0.9);
1749         hd0SPDany->Draw("sames");
1750         ctrsel->Update();
1751         TPaveStats *st2=(TPaveStats*)hd0SPDany->GetListOfFunctions()->FindObject("stats");
1752         st2->SetY1NDC(0.51);
1753         st2->SetY2NDC(0.7);
1754         st2->SetTextColor(4);
1755         hd0fb4->Draw("sames");
1756         ctrsel->Update();
1757         TPaveStats *st3=(TPaveStats*)hd0fb4->GetListOfFunctions()->FindObject("stats");
1758         st3->SetY1NDC(0.31);
1759         st3->SetY2NDC(0.5);
1760         st3->SetTextColor(2);
1761         hd0TPCITScuts->Draw("sames");
1762         ctrsel->Update();
1763         TPaveStats *st4=(TPaveStats*)hd0TPCITScuts->GetListOfFunctions()->FindObject("stats");
1764        st4->SetY1NDC(0.71);
1765        st4->SetY2NDC(0.9);
1766        st4->SetX1NDC(0.55);
1767        st4->SetX2NDC(0.75);
1768        st4->SetTextColor(kGreen+1);
1769        ctrsel->Modified();
1770        leg3->AddEntry(hd0SPD1,"kITSrefit+SPD inner","L");
1771        leg3->AddEntry(hd0SPDany,"kITSrefit+SPD any","L");
1772        leg3->AddEntry(hd0TPCITScuts,"TPC+ITS cuts+SPD any","L");
1773        leg3->AddEntry(hd0fb4,"Filter Bit 4","L");
1774        leg3->AddEntry(hd0SPD1, legend[i], "L");
1775       }
1776      else{
1777        hd0SPD1->SetStats(0);
1778        hd0SPD1->SetLineStyle(i+1);
1779        hd0SPD1->Draw("sames");
1780        hd0SPDany->SetStats(0);
1781        hd0SPDany->SetLineStyle(i+1);
1782        hd0TPCITScuts->SetStats(0);
1783        hd0TPCITScuts->SetLineStyle(i+1);
1784        hd0fb4->SetStats(0);
1785        hd0fb4->SetLineStyle(i+1);
1786        ctrsel->cd();
1787        hd0SPD1->Draw("sames");
1788        hd0SPDany->Draw("sames");
1789        hd0TPCITScuts->Draw("sames");
1790        hd0fb4->Draw("sames");
1791        leg3->AddEntry(hd0SPD1, legend[i], "L");
1792      }
1793
1794     }
1795   }
1796   leg3->Draw();
1797   ctrsel->SaveAs("ImpactParameterTrackSel.eps");
1798   ctrsel->SaveAs("ImpactParameterTrackSel.pdf");
1799 }
1800