PWG4/PartCorrDep/AliAnaCalorimeterQA.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15 /* $Id: $ */
  16
  17 //_________________________________________________________________________
  18 // Class to check results from simulations or reconstructed real data.
  19 // Fill few histograms and do some checking plots
  20 //
  21 //-- Author: Gustavo Conesa (INFN-LNF)
  22 //_________________________________________________________________________
  23
  24
  25 // --- ROOT system ---
  26 #include "Riostream.h"
  27 #include "TRefArray.h"
  28 #include "TParticle.h"
  29 #include "TCanvas.h"
  30 #include "TPad.h"
  31 #include "TROOT.h"
  32 //#include "TH3F.h"
  33 #include "TH2F.h"
  34 #include "TLegend.h"
  35
  36 //---- AliRoot system ----
  37 #include "AliAnaCalorimeterQA.h"
  38 #include "AliCaloTrackReader.h"
  39 #include "AliStack.h"
  40 #include "AliAODCaloCells.h"
  41 #include "AliAODCaloCluster.h"
  42 #include "AliFidutialCut.h"
  43 #include "AliESDtrack.h"
  44 #include "AliESDCaloCluster.h"
  45 #include "AliESDEvent.h"
  46
  47 ClassImp(AliAnaCalorimeterQA)
  48
  49 //____________________________________________________________________________
  50   AliAnaCalorimeterQA::AliAnaCalorimeterQA() :
  51     AliAnaPartCorrBaseClass(), fCalorimeter(""), fStyleMacro(""), fhE(0),fhPt(0),fhPhi(0),fhEta(0), fhEtaPhi(0),
  52     fhECharged(0),fhPtCharged(0),fhPhiCharged(0),fhEtaCharged(0), fhEtaPhiCharged(0),
  53     fhEChargedNoOut(0),fhPtChargedNoOut(0),fhPhiChargedNoOut(0),fhEtaChargedNoOut(0), fhEtaPhiChargedNoOut(0),
  54     fhDeltaE(0), fhDeltaPt(0),fhDeltaPhi(0),fhDeltaEta(0), fhRatioE(0), fhRatioPt(0),fhRatioPhi(0),fhRatioEta(0),
  55     fh2E(0),fh2Pt(0),fh2Phi(0),fh2Eta(0), fhIM(0), fhAsym(0), fhNCellsPerCluster(0), fhNClusters(0), fhNCells(0), fhAmplitude(0),
  56     fhGenGamPt(0),fhGenGamEta(0),fhGenGamPhi(0),fhGenPi0Pt(0),fhGenPi0Eta(0),fhGenPi0Phi(0),
  57     fhGenEtaPt(0),fhGenEtaEta(0),fhGenEtaPhi(0),fhGenOmegaPt(0),fhGenOmegaEta(0),fhGenOmegaPhi(0),
  58     fhGenElePt(0),fhGenEleEta(0),fhGenElePhi(0), fhEMVxyz(0),  fhEMR(0), fhHaVxyz(0),  fhHaR(0),
  59     fhGamE(0),fhGamPt(0),fhGamPhi(0),fhGamEta(0),
  60     fhGamDeltaE(0), fhGamDeltaPt(0),fhGamDeltaPhi(0),fhGamDeltaEta(0), fhGamRatioE(0), fhGamRatioPt(0),fhGamRatioPhi(0),fhGamRatioEta(0),
  61     fhGam2E(0),fhGam2Pt(0),fhGam2Phi(0),fhGam2Eta(0),
  62     fhEleE(0),fhElePt(0),fhElePhi(0),fhEleEta(0),
  63     fhPi0E(0),fhPi0Pt(0),fhPi0Phi(0),fhPi0Eta(0), fhNeHadE(0),fhNeHadPt(0),fhNeHadPhi(0),fhNeHadEta(0),
  64     fhChHadE(0),fhChHadPt(0),fhChHadPhi(0),fhChHadEta(0),
  65     fhGenGamAccE(0),fhGenGamAccPt(0),fhGenGamAccEta(0),fhGenGamAccPhi(0),
  66     fhGenPi0AccE(0),fhGenPi0AccPt(0),fhGenPi0AccEta(0),fhGenPi0AccPhi(0),
  67     fh1pOverE(0),fh1dR(0),fh2EledEdx(0), fh2MatchdEdx(0)
  68 {
  69   //Default Ctor
  70
  71   //Initialize parameters
  72   InitParameters();
  73 }
  74
  75 //____________________________________________________________________________
  76 AliAnaCalorimeterQA::AliAnaCalorimeterQA(const AliAnaCalorimeterQA & qa) :
  77   AliAnaPartCorrBaseClass(qa), fCalorimeter(qa.fCalorimeter), fStyleMacro(qa.fStyleMacro),
  78   fhE(qa.fhE),fhPt(qa.fhPt), fhPhi(qa.fhPhi), fhEta(qa.fhEta),
  79   fhEtaPhi(qa.fhEtaPhi), fhECharged(qa.fhECharged),fhPtCharged(qa.fhPtCharged),fhPhiCharged(qa.fhPhiCharged),
  80   fhEtaCharged(qa.fhEtaCharged), fhEtaPhiCharged(qa.fhEtaPhiCharged),
  81   fhEChargedNoOut(qa.fhEChargedNoOut),fhPtChargedNoOut(qa.fhPtChargedNoOut),fhPhiChargedNoOut(qa.fhPhiChargedNoOut),
  82   fhEtaChargedNoOut(qa.fhEtaChargedNoOut), fhEtaPhiChargedNoOut(qa.fhEtaPhiChargedNoOut),
  83   fhDeltaE(qa.fhDeltaE), fhDeltaPt(qa.fhDeltaPt), fhDeltaPhi(qa.fhDeltaPhi), fhDeltaEta(qa.fhDeltaEta),
  84   fhRatioE(qa.fhRatioE), fhRatioPt(qa.fhRatioPt), fhRatioPhi(qa.fhRatioPhi), fhRatioEta(qa.fhRatioEta),
  85   fh2E(qa.fh2E), fh2Pt(qa.fh2Pt), fh2Phi(qa.fh2Phi),fh2Eta(qa.fh2Eta),
  86   fhIM(qa.fhIM), fhAsym(qa.fhAsym), fhNCellsPerCluster(qa.fhNCellsPerCluster), fhNClusters(qa.fhNClusters),
  87   fhNCells(qa.fhNCells), fhAmplitude(qa.fhAmplitude),
  88   fhGenGamPt(qa.fhGenGamPt), fhGenGamEta(qa.fhGenGamEta), fhGenGamPhi(qa.fhGenGamPhi),
  89   fhGenPi0Pt(qa.fhGenPi0Pt), fhGenPi0Eta(qa.fhGenPi0Eta), fhGenPi0Phi(qa.fhGenPi0Phi),
  90   fhGenEtaPt(qa.fhGenEtaPt), fhGenEtaEta(qa.fhGenEtaEta), fhGenEtaPhi(qa.fhGenEtaPhi),
  91   fhGenOmegaPt(qa.fhGenOmegaPt), fhGenOmegaEta(qa.fhGenOmegaEta), fhGenOmegaPhi(qa.fhGenOmegaPhi),
  92   fhGenElePt(qa.fhGenElePt), fhGenEleEta(qa.fhGenEleEta), fhGenElePhi(qa.fhGenElePhi),
  93   fhEMVxyz(qa.fhEMVxyz),  fhEMR(qa.fhEMR), fhHaVxyz(qa.fhHaVxyz),  fhHaR(qa.fhHaR),
  94   fhGamE(qa.fhGamE),fhGamPt(qa.fhGamPt),fhGamPhi(qa.fhGamPhi),fhGamEta(qa.fhGamEta),
  95   fhGamDeltaE(qa.fhGamDeltaE), fhGamDeltaPt(qa.fhGamDeltaPt), fhGamDeltaPhi(qa.fhGamDeltaPhi), fhGamDeltaEta(qa.fhGamDeltaEta),
  96   fhGamRatioE(qa.fhGamRatioE), fhGamRatioPt(qa.fhGamRatioPt), fhGamRatioPhi(qa.fhGamRatioPhi), fhGamRatioEta(qa.fhGamRatioEta),
  97   fhGam2E(qa.fhGam2E), fhGam2Pt(qa.fhGam2Pt), fhGam2Phi(qa.fhGam2Phi),fhGam2Eta(qa.fhGam2Eta),
  98   fhEleE(qa.fhEleE),fhElePt(qa.fhElePt),fhElePhi(qa.fhElePhi),fhEleEta(qa.fhEleEta),
  99   fhPi0E(qa.fhPi0E),fhPi0Pt(qa.fhPi0Pt),fhPi0Phi(qa.fhPi0Phi),fhPi0Eta(qa.fhPi0Eta),
 100   fhNeHadE(qa.fhNeHadE),fhNeHadPt(qa.fhNeHadPt),fhNeHadPhi(qa.fhNeHadPhi),fhNeHadEta(qa.fhNeHadEta),
 101   fhChHadE(qa.fhChHadE),fhChHadPt(qa.fhChHadPt),fhChHadPhi(qa.fhChHadPhi),fhChHadEta(qa.fhChHadEta),
 102   fhGenGamAccE(qa.fhGenGamAccE),fhGenGamAccPt(qa.fhGenGamAccPt),fhGenGamAccEta(qa.fhGenGamAccEta),fhGenGamAccPhi(qa.fhGenGamAccPhi),
 103   fhGenPi0AccE(qa.fhGenPi0AccE),fhGenPi0AccPt(qa.fhGenPi0AccPt),fhGenPi0AccEta(qa.fhGenPi0AccEta),fhGenPi0AccPhi(qa.fhGenPi0AccPhi),
 104   fh1pOverE(qa.fh1pOverE),fh1dR(qa.fh1dR),fh2EledEdx(qa.fh2EledEdx), fh2MatchdEdx(qa.fh2MatchdEdx)
 105 {
 106   // cpy ctor
 107
 108 }
 109
 110 //_________________________________________________________________________
 111 AliAnaCalorimeterQA & AliAnaCalorimeterQA::operator = (const AliAnaCalorimeterQA & qa)
 112 {
 113   // assignment operator
 114
 115   if(this == &qa)return *this;
 116   ((AliAnaPartCorrBaseClass *)this)->operator=(qa);
 117
 118   fCalorimeter  = qa.fCalorimeter;
 119   fStyleMacro   = qa.fStyleMacro;
 120
 121   fhE      = qa.fhE;
 122   fhPt     = qa.fhPt;
 123   fhPhi    = qa.fhPhi;
 124   fhEta    = qa.fhEta;
 125   fhEtaPhi = qa.fhEtaPhi;
 126
 127   fhECharged      = qa.fhECharged;
 128   fhPtCharged     = qa.fhPtCharged;
 129   fhPhiCharged    = qa.fhPhiCharged;
 130   fhEtaCharged    = qa.fhEtaCharged;
 131   fhEtaPhiCharged = qa.fhEtaPhiCharged;
 132
 133   fhEChargedNoOut      = qa.fhEChargedNoOut;
 134   fhPtChargedNoOut     = qa.fhPtChargedNoOut;
 135   fhPhiChargedNoOut    = qa.fhPhiChargedNoOut;
 136   fhEtaChargedNoOut    = qa.fhEtaChargedNoOut;
 137   fhEtaPhiChargedNoOut = qa.fhEtaPhiChargedNoOut;
 138
 139   fhIM   = qa.fhIM;
 140   fhAsym = qa.fhAsym;
 141
 142   fhNCellsPerCluster = qa.fhNCellsPerCluster;
 143   fhNClusters        = qa.fhNClusters;
 144
 145   fhDeltaE   = qa.fhDeltaE;
 146   fhDeltaPt  = qa.fhDeltaPt;
 147   fhDeltaPhi = qa.fhDeltaPhi;
 148   fhDeltaEta = qa.fhDeltaEta;
 149
 150   fhRatioE   = qa.fhRatioE;
 151   fhRatioPt  = qa.fhRatioPt;
 152   fhRatioPhi = qa.fhRatioPhi;
 153   fhRatioEta = qa.fhRatioEta;
 154
 155
 156   fh2E   = qa.fh2E;
 157   fh2Pt  = qa.fh2Pt;
 158   fh2Phi = qa.fh2Phi;
 159   fh2Eta = qa.fh2Eta;
 160
 161   fhNCells    = qa.fhNCells;
 162   fhAmplitude = qa.fhAmplitude;
 163
 164   fhGenGamPt = qa.fhGenGamPt  ; fhGenGamEta = qa.fhGenGamEta  ; fhGenGamPhi = qa.fhGenGamPhi  ;
 165   fhGenPi0Pt = qa.fhGenPi0Pt  ; fhGenPi0Eta = qa.fhGenPi0Eta  ; fhGenPi0Phi = qa.fhGenPi0Phi  ;
 166   fhGenEtaPt = qa.fhGenEtaPt  ; fhGenEtaEta = qa.fhGenEtaEta  ; fhGenEtaPhi = qa.fhGenEtaPhi  ;
 167   fhGenOmegaPt = qa.fhGenOmegaPt  ; fhGenOmegaEta = qa.fhGenOmegaEta  ; fhGenOmegaPhi = qa.fhGenOmegaPhi  ;
 168   fhGenElePt = qa.fhGenElePt  ; fhGenEleEta = qa.fhGenEleEta  ; fhGenElePhi = qa.fhGenElePhi ;
 169
 170   fhEMVxyz = qa.fhEMVxyz ;  fhEMR = qa.fhEMR ; fhHaVxyz = qa.fhHaVxyz ;  fhHaR = qa.fhHaR ;
 171   fhGamE = qa.fhGamE ;fhGamPt = qa.fhGamPt ;fhGamPhi = qa.fhGamPhi ;fhGamEta = qa.fhGamEta ;
 172   fhGamDeltaE   = qa.fhDeltaE; fhGamDeltaPt  = qa.fhDeltaPt; fhGamDeltaPhi = qa.fhDeltaPhi; fhGamDeltaEta = qa.fhDeltaEta;
 173
 174   fhGamRatioE   = qa.fhGamRatioE;  fhGamRatioPt  = qa.fhGamRatioPt;  fhGamRatioPhi = qa.fhGamRatioPhi;  fhGamRatioEta = qa.fhGamRatioEta;
 175   fhGam2E   = qa.fhGam2E;        fhGam2Pt  = qa.fhGam2Pt; fhGam2Phi = qa.fhGam2Phi; fhGam2Eta = qa.fhGam2Eta;
 176
 177   fhEleE = qa.fhEleE ;fhElePt = qa.fhElePt ;fhElePhi = qa.fhElePhi ;fhEleEta = qa.fhEleEta ;
 178   fhPi0E = qa.fhPi0E ;fhPi0Pt = qa.fhPi0Pt ;fhPi0Phi = qa.fhPi0Phi ;fhPi0Eta = qa.fhPi0Eta ;
 179   fhNeHadE = qa.fhNeHadE ;fhNeHadPt = qa.fhNeHadPt ;fhNeHadPhi = qa.fhNeHadPhi ;fhNeHadEta = qa.fhNeHadEta ;
 180   fhChHadE = qa.fhChHadE ;fhChHadPt = qa.fhChHadPt ;fhChHadPhi = qa.fhChHadPhi ;fhChHadEta = qa.fhChHadEta ;
 181   fhGenGamAccE = qa.fhGenGamAccE ;  fhGenPi0AccE = qa.fhGenPi0AccE ;
 182   fhGenGamAccPt = qa.fhGenGamAccPt ;fhGenGamAccEta = qa.fhGenGamAccEta ;fhGenGamAccPhi = qa.fhGenGamAccPhi ;
 183   fhGenPi0AccPt = qa.fhGenPi0AccPt ;fhGenPi0AccEta = qa.fhGenPi0AccEta; fhGenPi0AccPhi = qa.fhGenPi0AccPhi ;
 184
 185   fh1pOverE = qa.fh1pOverE;
 186   fh1dR = qa.fh1dR;
 187   fh2MatchdEdx = qa.fh2MatchdEdx;
 188   fh2EledEdx = qa.fh2EledEdx;
 189
 190   return *this;
 191
 192 }
 193
 194 //________________________________________________________________________
 195 TList *  AliAnaCalorimeterQA::GetCreateOutputObjects()
 196 {
 197         // Create histograms to be saved in output file and
 198         // store them in fOutputContainer
 199
 200         TList * outputContainer = new TList() ;
 201         outputContainer->SetName("ExampleHistos") ;
 202
 203         Int_t nptbins  = GetHistoNPtBins();
 204         Int_t nphibins = GetHistoNPhiBins();
 205         Int_t netabins = GetHistoNEtaBins();
 206         Float_t ptmax  = GetHistoPtMax();
 207         Float_t phimax = GetHistoPhiMax();
 208         Float_t etamax = GetHistoEtaMax();
 209         Float_t ptmin  = GetHistoPtMin();
 210         Float_t phimin = GetHistoPhiMin();
 211         Float_t etamin = GetHistoEtaMin();
 212
 213         fhE  = new TH1F ("hE","E reconstructed", nptbins,ptmin,ptmax);
 214         fhE->SetXTitle("E (GeV)");
 215         outputContainer->Add(fhE);
 216
 217         fhPt  = new TH1F ("hPt","p_{T} reconstructed", nptbins,ptmin,ptmax);
 218         fhPt->SetXTitle("p_{T} (GeV/c)");
 219         outputContainer->Add(fhPt);
 220
 221         fhPhi  = new TH1F ("hPhi","#phi reconstructed",nphibins,phimin,phimax);
 222         fhPhi->SetXTitle("#phi (rad)");
 223         outputContainer->Add(fhPhi);
 224
 225         fhEta  = new TH1F ("hEta","#eta reconstructed",netabins,etamin,etamax);
 226         fhEta->SetXTitle("#eta ");
 227         outputContainer->Add(fhEta);
 228
 229         fhEtaPhi  = new TH2F ("hEtaPhi","#eta vs #phi, reconstructed",netabins,etamin,etamax,nphibins,phimin,phimax);
 230         fhEtaPhi->SetXTitle("#eta ");
 231         fhEtaPhi->SetYTitle("#phi ");
 232         outputContainer->Add(fhEtaPhi);
 233
 234         fhECharged  = new TH1F ("hECharged","E reconstructed, matched with track", nptbins,ptmin,ptmax);
 235         fhECharged->SetXTitle("E (GeV)");
 236         outputContainer->Add(fhECharged);
 237
 238         fhPtCharged  = new TH1F ("hPtCharged","p_{T} reconstructed, matched with track", nptbins,ptmin,ptmax);
 239         fhPtCharged->SetXTitle("p_{T} (GeV/c)");
 240         outputContainer->Add(fhPtCharged);
 241
 242         fhPhiCharged  = new TH1F ("hPhiCharged","#phi reconstructed, matched with track",nphibins,phimin,phimax);
 243         fhPhiCharged->SetXTitle("#phi (rad)");
 244         outputContainer->Add(fhPhiCharged);
 245
 246         fhEtaCharged  = new TH1F ("hEtaCharged","#eta reconstructed, matched with track",netabins,etamin,etamax);
 247         fhEtaCharged->SetXTitle("#eta ");
 248         outputContainer->Add(fhEtaCharged);
 249
 250         fhEtaPhiCharged  = new TH2F ("hEtaPhiCharged","#eta vs #phi, reconstructed , matched with track",netabins,etamin,etamax,nphibins,phimin,phimax);
 251         fhEtaPhiCharged->SetXTitle("#eta ");
 252         fhEtaPhiCharged->SetYTitle("#phi ");
 253         outputContainer->Add(fhEtaPhiCharged);
 254
 255
 256         fhEChargedNoOut  = new TH1F ("hEChargedNoOut","E reconstructed, matched with track, no output track params", nptbins,ptmin,ptmax);
 257         fhEChargedNoOut->SetXTitle("E (GeV)");
 258         outputContainer->Add(fhEChargedNoOut);
 259
 260         fhPtChargedNoOut  = new TH1F ("hPtChargedNoOut","p_{T} reconstructed, matched with track, no output track params", nptbins,ptmin,ptmax);
 261         fhPtChargedNoOut->SetXTitle("p_{T} (GeV/c)");
 262         outputContainer->Add(fhPtChargedNoOut);
 263
 264         fhPhiChargedNoOut  = new TH1F ("hPhiChargedNoOut","#phi reconstructed, matched with track, no output track params",nphibins,phimin,phimax);
 265         fhPhiChargedNoOut->SetXTitle("#phi (rad)");
 266         outputContainer->Add(fhPhiChargedNoOut);
 267
 268         fhEtaChargedNoOut  = new TH1F ("hEtaChargedNoOut","#eta reconstructed, matched with track, no output track params",netabins,etamin,etamax);
 269         fhEtaChargedNoOut->SetXTitle("#eta ");
 270         outputContainer->Add(fhEtaChargedNoOut);
 271
 272         fhEtaPhiChargedNoOut  = new TH2F ("hEtaPhiChargedNoOut","#eta vs #phi, reconstructed , matched with track, no output track params",netabins,etamin,etamax,nphibins,phimin,phimax);
 273         fhEtaPhiChargedNoOut->SetXTitle("#eta ");
 274         fhEtaPhiChargedNoOut->SetYTitle("#phi ");
 275         outputContainer->Add(fhEtaPhiChargedNoOut);
 276
 277         fhIM  = new TH2F ("hIM","Cluster pairs Invariant mass vs reconstructed pair energy",nptbins,ptmin,ptmax,200,0,1);
 278         fhIM->SetXTitle("E_{cluster pairs} (GeV) ");
 279         fhIM->SetYTitle("M_{cluster pairs} (GeV/c^{2})");
 280         outputContainer->Add(fhIM);
 281
 282         fhAsym  = new TH2F ("hAssym","Cluster pairs Asymmetry vs reconstructed pair energy",nptbins,ptmin,ptmax,100,0,1);
 283         fhAsym->SetXTitle("E_{cluster pairs} (GeV) ");
 284         fhAsym->SetYTitle("Asymmetry");
 285         outputContainer->Add(fhAsym);
 286
 287         fhNCellsPerCluster  = new TH2F ("hNCellsPerCluster","# cells per cluster vs cluster energy", nptbins,ptmin,ptmax, 100,0,1000);
 288         fhNCellsPerCluster->SetXTitle("E (GeV)");
 289         fhNCellsPerCluster->SetYTitle("n cells");
 290         outputContainer->Add(fhNCellsPerCluster);
 291
 292         fhNClusters  = new TH1F ("hNClusters","# clusters", 300,0,300);
 293         fhNClusters->SetXTitle("number of clusters");
 294         outputContainer->Add(fhNClusters);
 295
 296         //Calo cells
 297         fhNCells  = new TH1F ("hNCells","# cells", 13000,0,13000);
 298         fhNCells->SetXTitle("n cells");
 299         outputContainer->Add(fhNCells);
 300
 301         fhAmplitude  = new TH1F ("hAmplitude","Amplitude", 100,0,1000);
 302         fhAmplitude->SetXTitle("Amplitude ");
 303         outputContainer->Add(fhAmplitude);
 304
 305         if(IsDataMC()){
 306
 307                 fhDeltaE  = new TH1F ("hDeltaE","MC - Reco E ", 200,-50,50);
 308                 fhDeltaE->SetXTitle("#Delta E (GeV)");
 309                 outputContainer->Add(fhDeltaE);
 310
 311                 fhDeltaPt  = new TH1F ("hDeltaPt","MC - Reco p_{T} ", 200,-50,50);
 312                 fhDeltaPt->SetXTitle("#Delta p_{T} (GeV/c)");
 313                 outputContainer->Add(fhDeltaPt);
 314
 315                 fhDeltaPhi  = new TH1F ("hDeltaPhi","MC - Reco #phi ",100,-2,2);
 316                 fhDeltaPhi->SetXTitle("#Delta #phi (rad)");
 317                 outputContainer->Add(fhDeltaPhi);
 318
 319                 fhDeltaEta  = new TH1F ("hDeltaEta","MC- Reco #eta",100,-1,1);
 320                 fhDeltaEta->SetXTitle("#Delta #eta ");
 321                 outputContainer->Add(fhDeltaEta);
 322
 323                 fhRatioE  = new TH1F ("hRatioE","Reco/MC E ", 200,0,2);
 324                 fhRatioE->SetXTitle("E_{reco}/E_{gen}");
 325                 outputContainer->Add(fhRatioE);
 326
 327                 fhRatioPt  = new TH1F ("hRatioPt","Reco/MC p_{T} ", 200,0,2);
 328                 fhRatioPt->SetXTitle("p_{T, reco}/p_{T, gen}");
 329                 outputContainer->Add(fhRatioPt);
 330
 331                 fhRatioPhi  = new TH1F ("hRatioPhi","Reco/MC #phi ",200,0,2);
 332                 fhRatioPhi->SetXTitle("#phi_{reco}/#phi_{gen}");
 333                 outputContainer->Add(fhRatioPhi);
 334
 335                 fhRatioEta  = new TH1F ("hRatioEta","Reco/MC #eta",200,0,2);
 336                 fhRatioEta->SetXTitle("#eta_{reco}/#eta_{gen} ");
 337                 outputContainer->Add(fhRatioEta);
 338
 339                 fh2E  = new TH2F ("h2E","E distribution, reconstructed vs generated", nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
 340                 fh2E->SetXTitle("E_{rec} (GeV)");
 341                 fh2E->SetYTitle("E_{gen} (GeV)");
 342                 outputContainer->Add(fh2E);
 343
 344                 fh2Pt  = new TH2F ("h2Pt","p_T distribution, reconstructed vs generated", nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
 345                 fh2Pt->SetXTitle("p_{T,rec} (GeV/c)");
 346                 fh2Pt->SetYTitle("p_{T,gen} (GeV/c)");
 347                 outputContainer->Add(fh2Pt);
 348
 349                 fh2Phi  = new TH2F ("h2Phi","#phi distribution, reconstructed vs generated", nphibins,phimin,phimax, nphibins,phimin,phimax);
 350                 fh2Phi->SetXTitle("#phi_{rec} (rad)");
 351                 fh2Phi->SetYTitle("#phi_{gen} (rad)");
 352                 outputContainer->Add(fh2Phi);
 353
 354                 fh2Eta  = new TH2F ("h2Eta","#eta distribution, reconstructed vs generated", netabins,etamin,etamax,netabins,etamin,etamax);
 355                 fh2Eta->SetXTitle("#eta_{rec} ");
 356                 fh2Eta->SetYTitle("#eta_{gen} ");
 357                 outputContainer->Add(fh2Eta);
 358
 359                 //Fill histos depending on origin of cluster
 360                 fhGamE  = new TH2F ("hGamE","E reconstructed vs E generated from #gamma", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
 361                 fhGamE->SetXTitle("E (GeV)");
 362                 outputContainer->Add(fhGamE);
 363
 364                 fhGamPt  = new TH2F ("hGamPt","p_{T} reconstructed vs E generated from #gamma", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
 365                 fhGamPt->SetXTitle("p_{T} (GeV/c)");
 366                 outputContainer->Add(fhGamPt);
 367
 368                 fhGamPhi  = new TH2F ("hGamPhi","#phi reconstructed vs E generated from #gamma",nphibins,phimin,phimax,nphibins,phimin,phimax);
 369                 fhGamPhi->SetXTitle("#phi (rad)");
 370                 outputContainer->Add(fhGamPhi);
 371
 372                 fhGamEta  = new TH2F ("hGamEta","#eta reconstructed vs E generated from #gamma",netabins,etamin,etamax,netabins,etamin,etamax);
 373                 fhGamEta->SetXTitle("#eta ");
 374                 outputContainer->Add(fhGamEta);
 375
 376                 fhGamDeltaE  = new TH1F ("hGamDeltaE","#gamma MC - Reco E ", 200,-50,50);
 377                 fhGamDeltaE->SetXTitle("#Delta E (GeV)");
 378                 outputContainer->Add(fhGamDeltaE);
 379
 380                 fhGamDeltaPt  = new TH1F ("hGamDeltaPt","#gamma MC - Reco p_{T} ", 200,-50,50);
 381                 fhGamDeltaPt->SetXTitle("#Delta p_{T} (GeV/c)");
 382                 outputContainer->Add(fhGamDeltaPt);
 383
 384                 fhGamDeltaPhi  = new TH1F ("hGamDeltaPhi","#gamma MC - Reco #phi ",100,-2,2);
 385                 fhGamDeltaPhi->SetXTitle("#Delta #phi (rad)");
 386                 outputContainer->Add(fhGamDeltaPhi);
 387
 388                 fhGamDeltaEta  = new TH1F ("hGamDeltaEta","#gamma MC- Reco #eta",100,-1,1);
 389                 fhGamDeltaEta->SetXTitle("#Delta #eta ");
 390                 outputContainer->Add(fhGamDeltaEta);
 391
 392                 fhGamRatioE  = new TH1F ("hGamRatioE","#gamma Reco/MC E ", 200,0,2);
 393                 fhGamRatioE->SetXTitle("E_{reco}/E_{gen}");
 394                 outputContainer->Add(fhGamRatioE);
 395
 396                 fhGamRatioPt  = new TH1F ("hGamRatioPt","#gamma Reco/MC p_{T} ", 200,0,2);
 397                 fhGamRatioPt->SetXTitle("p_{T, reco}/p_{T, gen}");
 398                 outputContainer->Add(fhGamRatioPt);
 399
 400                 fhGamRatioPhi  = new TH1F ("hGamRatioPhi","#gamma Reco/MC #phi ",200,0,2);
 401                 fhGamRatioPhi->SetXTitle("#phi_{reco}/#phi_{gen}");
 402                 outputContainer->Add(fhGamRatioPhi);
 403
 404                 fhGamRatioEta  = new TH1F ("hGamRatioEta","#gamma Reco/MC #eta",200,0,2);
 405                 fhGamRatioEta->SetXTitle("#eta_{reco}/#eta_{gen} ");
 406                 outputContainer->Add(fhGamRatioEta);
 407
 408                 fhGam2E  = new TH2F ("hGam2E","#gamma E distribution, reconstructed vs generated", nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
 409                 fhGam2E->SetXTitle("E_{rec} (GeV)");
 410                 fhGam2E->SetYTitle("E_{gen} (GeV)");
 411                 outputContainer->Add(fhGam2E);
 412
 413                 fhGam2Pt  = new TH2F ("hGam2Pt","#gamma p_T distribution, reconstructed vs generated", nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
 414                 fhGam2Pt->SetXTitle("p_{T,rec} (GeV/c)");
 415                 fhGam2Pt->SetYTitle("p_{T,gen} (GeV/c)");
 416                 outputContainer->Add(fhGam2Pt);
 417
 418                 fhGam2Phi  = new TH2F ("hGam2Phi","#gamma #phi distribution, reconstructed vs generated", nphibins,phimin,phimax, nphibins,phimin,phimax);
 419                 fhGam2Phi->SetXTitle("#phi_{rec} (rad)");
 420                 fhGam2Phi->SetYTitle("#phi_{gen} (rad)");
 421                 outputContainer->Add(fhGam2Phi);
 422
 423                 fhGam2Eta  = new TH2F ("hGam2Eta","#gamma #eta distribution, reconstructed vs generated", netabins,etamin,etamax,netabins,etamin,etamax);
 424                 fhGam2Eta->SetXTitle("#eta_{rec} ");
 425                 fhGam2Eta->SetYTitle("#eta_{gen} ");
 426                 outputContainer->Add(fhGam2Eta);
 427
 428                 fhPi0E  = new TH2F ("hPi0E","E reconstructed vs E generated from #pi^0", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
 429                 fhPi0E->SetXTitle("E (GeV)");
 430                 outputContainer->Add(fhPi0E);
 431
 432                 fhPi0Pt  = new TH2F ("hPi0Pt","p_{T} reconstructed vs E generated from #pi^{0}", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
 433                 fhPi0Pt->SetXTitle("p_{T} (GeV/c)");
 434                 outputContainer->Add(fhPi0Pt);
 435
 436                 fhPi0Phi  = new TH2F ("hPi0Phi","#phi reconstructed vs E generated from #pi^{0}",nphibins,phimin,phimax,nphibins,phimin,phimax);
 437                 fhPi0Phi->SetXTitle("#phi (rad)");
 438                 outputContainer->Add(fhPi0Phi);
 439
 440                 fhPi0Eta  = new TH2F ("hPi0Eta","#eta reconstructed vs E generated from #pi^{0}",netabins,etamin,etamax,netabins,etamin,etamax);
 441                 fhPi0Eta->SetXTitle("#eta ");
 442                 outputContainer->Add(fhPi0Eta);
 443
 444                 fhEleE  = new TH2F ("hEleE","E reconstructed vs E generated from e^{#pm}", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
 445                 fhEleE->SetXTitle("E (GeV)");
 446                 outputContainer->Add(fhEleE);
 447
 448                 fhElePt  = new TH2F ("hElePt","p_{T} reconstructed vs E generated from e^{#pm}", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
 449                 fhElePt->SetXTitle("p_{T} (GeV/c)");
 450                 outputContainer->Add(fhElePt);
 451
 452                 fhElePhi  = new TH2F ("hElePhi","#phi reconstructed vs E generated frome^{#pm}",nphibins,phimin,phimax,nphibins,phimin,phimax);
 453                 fhElePhi->SetXTitle("#phi (rad)");
 454                 outputContainer->Add(fhElePhi);
 455
 456                 fhEleEta  = new TH2F ("hEleEta","#eta reconstructed vs E generated from e^{#pm}",netabins,etamin,etamax,netabins,etamin,etamax);
 457                 fhEleEta->SetXTitle("#eta ");
 458                 outputContainer->Add(fhEleEta);
 459
 460                 fhNeHadE  = new TH2F ("hNeHadE","E reconstructed vs E generated from #pi^0", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
 461                 fhNeHadE->SetXTitle("E (GeV)");
 462                 outputContainer->Add(fhNeHadE);
 463
 464                 fhNeHadPt  = new TH2F ("hNeHadPt","p_{T} reconstructed vs E generated from neutral hadron", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
 465                 fhNeHadPt->SetXTitle("p_{T} (GeV/c)");
 466                 outputContainer->Add(fhNeHadPt);
 467
 468                 fhNeHadPhi  = new TH2F ("hNeHadPhi","#phi reconstructed vs E generated from neutral hadron",nphibins,phimin,phimax,nphibins,phimin,phimax);
 469                 fhNeHadPhi->SetXTitle("#phi (rad)");
 470                 outputContainer->Add(fhNeHadPhi);
 471
 472                 fhNeHadEta  = new TH2F ("hNeHadEta","#eta reconstructed vs E generated from neutral hadron",netabins,etamin,etamax,netabins,etamin,etamax);
 473                 fhNeHadEta->SetXTitle("#eta ");
 474                 outputContainer->Add(fhNeHadEta);
 475
 476                 fhChHadE  = new TH2F ("hChHadE","E reconstructed vs E generated from #pi^0", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
 477                 fhChHadE->SetXTitle("E (GeV)");
 478                 outputContainer->Add(fhChHadE);
 479
 480                 fhChHadPt  = new TH2F ("hChHadPt","p_{T} reconstructed vs E generated from charged hadron", nptbins,ptmin,ptmax, nptbins,ptmin,ptmax);
 481                 fhChHadPt->SetXTitle("p_{T} (GeV/c)");
 482                 outputContainer->Add(fhChHadPt);
 483
 484                 fhChHadPhi  = new TH2F ("hChHadPhi","#phi reconstructed vs E generated from charged hadron",nphibins,phimin,phimax,nphibins,phimin,phimax);
 485                 fhChHadPhi->SetXTitle("#phi (rad)");
 486                 outputContainer->Add(fhChHadPhi);
 487
 488                 fhChHadEta  = new TH2F ("hChHadEta","#eta reconstructed vs E generated from charged hadron",netabins,etamin,etamax,netabins,etamin,etamax);
 489                 fhChHadEta->SetXTitle("#eta ");
 490                 outputContainer->Add(fhChHadEta);
 491
 492                 //Vertex of generated particles
 493
 494                 fhEMVxyz  = new TH2F ("hEMVxyz","Production vertex of reconstructed ElectroMagnetic particles",100,0,500,100,0,500);//,100,0,500);
 495                 fhEMVxyz->SetXTitle("v_{x}");
 496                 fhEMVxyz->SetYTitle("v_{y}");
 497                 //fhEMVxyz->SetZTitle("v_{z}");
 498                 outputContainer->Add(fhEMVxyz);
 499
 500                 fhHaVxyz  = new TH2F ("hHaVxyz","Production vertex of reconstructed hadrons",100,0,500,100,0,500);//,100,0,500);
 501                 fhHaVxyz->SetXTitle("v_{x}");
 502                 fhHaVxyz->SetYTitle("v_{y}");
 503                 //fhHaVxyz->SetZTitle("v_{z}");
 504                 outputContainer->Add(fhHaVxyz);
 505
 506                 fhEMR  = new TH2F ("hEMR","Distance to production vertex of reconstructed ElectroMagnetic particles vs E rec",nptbins,ptmin,ptmax,100,0,500);
 507                 fhEMR->SetXTitle("E (GeV)");
 508                 fhEMR->SetYTitle("TMath::Sqrt(v_{x}^{2}+v_{y}^{2})");
 509                 outputContainer->Add(fhEMR);
 510
 511                 fhHaR  = new TH2F ("hHaR","Distance to production vertex of reconstructed Hadrons vs E rec",nptbins,ptmin,ptmax,100,0,500);
 512                 fhHaR->SetXTitle("E (GeV)");
 513                 fhHaR->SetYTitle("TMath::Sqrt(v_{x}^{2}+v_{y}^{2})");
 514                 outputContainer->Add(fhHaR);
 515
 516
 517                 //Pure MC
 518                 fhGenGamPt  = new TH1F("hGenGamPt" ,"p_{T} of generated #gamma",nptbins,ptmin,ptmax);
 519                 fhGenGamEta = new TH1F("hGenGamEta","Y of generated #gamma",netabins,etamin,etamax);
 520                 fhGenGamPhi = new TH1F("hGenGamPhi","#phi of generated #gamma",nphibins,phimin,phimax);
 521
 522                 fhGenPi0Pt  = new TH1F("hGenPi0Pt" ,"p_{T} of generated #pi^{0}",nptbins,ptmin,ptmax);
 523                 fhGenPi0Eta = new TH1F("hGenPi0Eta","Y of generated #pi^{0}",netabins,etamin,etamax);
 524                 fhGenPi0Phi = new TH1F("hGenPi0Phi","#phi of generated #pi^{0}",nphibins,phimin,phimax);
 525
 526                 fhGenEtaPt  = new TH1F("hGenEtaPt" ,"p_{T} of generated #eta",nptbins,ptmin,ptmax);
 527                 fhGenEtaEta = new TH1F("hGenEtaEta","Y of generated #eta",netabins,etamin,etamax);
 528                 fhGenEtaPhi = new TH1F("hGenEtaPhi","#phi of generated #eta",nphibins,phimin,phimax);
 529
 530                 fhGenOmegaPt  = new TH1F("hGenOmegaPt" ,"p_{T} of generated #omega",nptbins,ptmin,ptmax);
 531                 fhGenOmegaEta = new TH1F("hGenOmegaEta","Y of generated #omega",netabins,etamin,etamax);
 532                 fhGenOmegaPhi = new TH1F("hGenOmegaPhi","#phi of generated #omega",nphibins,phimin,phimax);
 533
 534                 fhGenElePt  = new TH1F("hGenElePt" ,"p_{T} of generated e^{#pm}",nptbins,ptmin,ptmax);
 535                 fhGenEleEta = new TH1F("hGenEleEta","Y of generated  e^{#pm}",netabins,etamin,etamax);
 536                 fhGenElePhi = new TH1F("hGenElePhi","#phi of generated  e^{#pm}",nphibins,phimin,phimax);
 537
 538                 fhGenGamPt->SetXTitle("p_{T} (GeV/c)");
 539                 fhGenGamEta->SetXTitle("#eta");
 540                 fhGenGamPhi->SetXTitle("#phi (rad)");
 541                 outputContainer->Add(fhGenGamPt);
 542                 outputContainer->Add(fhGenGamEta);
 543                 outputContainer->Add(fhGenGamPhi);
 544
 545                 fhGenPi0Pt->SetXTitle("p_{T} (GeV/c)");
 546                 fhGenPi0Eta->SetXTitle("#eta");
 547                 fhGenPi0Phi->SetXTitle("#phi (rad)");
 548                 outputContainer->Add(fhGenPi0Pt);
 549                 outputContainer->Add(fhGenPi0Eta);
 550                 outputContainer->Add(fhGenPi0Phi);
 551
 552                 fhGenEtaPt->SetXTitle("p_{T} (GeV/c)");
 553                 fhGenEtaEta->SetXTitle("#eta");
 554                 fhGenEtaPhi->SetXTitle("#phi (rad)");
 555                 outputContainer->Add(fhGenEtaPt);
 556                 outputContainer->Add(fhGenEtaEta);
 557                 outputContainer->Add(fhGenEtaPhi);
 558
 559                 fhGenOmegaPt->SetXTitle("p_{T} (GeV/c)");
 560                 fhGenOmegaEta->SetXTitle("#eta");
 561                 fhGenOmegaPhi->SetXTitle("#phi (rad)");
 562                 outputContainer->Add(fhGenOmegaPt);
 563                 outputContainer->Add(fhGenOmegaEta);
 564                 outputContainer->Add(fhGenOmegaPhi);
 565
 566                 fhGenElePt->SetXTitle("p_{T} (GeV/c)");
 567                 fhGenEleEta->SetXTitle("#eta");
 568                 fhGenElePhi->SetXTitle("#phi (rad)");
 569                 outputContainer->Add(fhGenElePt);
 570                 outputContainer->Add(fhGenEleEta);
 571                 outputContainer->Add(fhGenElePhi);
 572
 573                 fhGenGamAccE   = new TH1F("hGenGamAccE" ,"E of generated #gamma in calorimeter acceptance",nptbins,ptmin,ptmax);
 574                 fhGenGamAccPt  = new TH1F("hGenGamAccPt" ,"p_{T} of generated #gamma in calorimeter acceptance",nptbins,ptmin,ptmax);
 575                 fhGenGamAccEta = new TH1F("hGenGamAccEta","Y of generated #gamma in calorimeter acceptance",netabins,etamin,etamax);
 576                 fhGenGamAccPhi = new TH1F("hGenGamAccPhi","#phi of generated #gamma  in calorimeter acceptance",nphibins,phimin,phimax);
 577
 578                 fhGenPi0AccE  = new TH1F("hGenPi0AccE" ,"E of generated #pi^{0} in calorimeter acceptance",nptbins,ptmin,ptmax);
 579                 fhGenPi0AccPt  = new TH1F("hGenPi0AccPt" ,"p_{T} of generated #pi^{0} in calorimeter acceptance",nptbins,ptmin,ptmax);
 580                 fhGenPi0AccEta = new TH1F("hGenPi0AccEta","Y of generated #pi^{0} in calorimeter acceptance",netabins,etamin,etamax);
 581                 fhGenPi0AccPhi = new TH1F("hGenPi0AccPhi","#phi of generated #pi^{0} in calorimeter acceptance",nphibins,phimin,phimax);
 582
 583                 fhGenGamAccE  ->SetXTitle("E (GeV)");
 584                 fhGenGamAccPt ->SetXTitle("p_{T} (GeV/c)");
 585                 fhGenGamAccEta->SetXTitle("#eta");
 586                 fhGenGamAccPhi->SetXTitle("#phi (rad)");
 587                 outputContainer->Add(fhGenGamAccE);
 588                 outputContainer->Add(fhGenGamAccPt);
 589                 outputContainer->Add(fhGenGamAccEta);
 590                 outputContainer->Add(fhGenGamAccPhi);
 591
 592                 fhGenPi0AccE  ->SetXTitle("E (GeV)");
 593                 fhGenPi0AccPt ->SetXTitle("p_{T} (GeV/c)");
 594                 fhGenPi0AccEta->SetXTitle("#eta");
 595                 fhGenPi0AccPhi->SetXTitle("#phi (rad)");
 596                 outputContainer->Add(fhGenPi0AccE);
 597                 outputContainer->Add(fhGenPi0AccPt);
 598                 outputContainer->Add(fhGenPi0AccEta);
 599                 outputContainer->Add(fhGenPi0AccPhi);
 600
 601         }
 602
 603
 604         fh1pOverE = new TH1F("h1pOverE","EMCAL-TRACK matches p/E",100,0.,10.);
 605         fh1pOverE->SetXTitle("p/E");
 606         outputContainer->Add(fh1pOverE);
 607
 608         fh1dR = new TH1F("h1dR","EMCAL-TRACK matches dR",300, 0.,TMath::Pi());
 609         fh1dR->SetXTitle("#Delta R (rad)");
 610         outputContainer->Add(fh1dR) ;
 611
 612         fh2MatchdEdx = new TH2F("h2MatchdEdx","dE/dx vs. p for all matches",200,0.,50.,200,0.,400.);
 613         fh2MatchdEdx->SetXTitle("p (GeV/c)");
 614         fh2MatchdEdx->SetYTitle("<dE/dx>");
 615         outputContainer->Add(fh2MatchdEdx);
 616
 617         fh2EledEdx = new TH2F("h2EledEdx","dE/dx vs. p for electrons",200,0.,50.,200,0.,400.);
 618         fh2EledEdx->SetXTitle("p (GeV/c)");
 619         fh2EledEdx->SetYTitle("<dE/dx>");
 620         outputContainer->Add(fh2EledEdx) ;
 621
 622         return outputContainer;
 623 }
 624
 625 //__________________________________________________
 626 void AliAnaCalorimeterQA::Init()
 627 {
 628         //Check if the data or settings are ok
 629         if(fCalorimeter != "PHOS" && fCalorimeter !="EMCAL"){
 630                 printf("AliAnaCalorimeterQA::Init() - Wrong calorimeter name <%s>, END\n", fCalorimeter.Data());
 631                 abort();
 632         }
 633
 634         if(GetReader()->GetDataType()== AliCaloTrackReader::kMC){
 635                 printf("AliAnaCalorimeterQA::Init() - Analysis of reconstructed data, MC reader not aplicable\n");
 636                 abort();
 637         }
 638
 639 }
 640
 641
 642 //__________________________________________________
 643 void AliAnaCalorimeterQA::InitParameters()
 644 {
 645   //Initialize the parameters of the analysis.
 646   AddToHistogramsName("AnaCaloQA_");
 647
 648   fCalorimeter = "EMCAL"; //or PHOS
 649   fStyleMacro = "" ;
 650 }
 651
 652 //__________________________________________________________________
 653 void AliAnaCalorimeterQA::Print(const Option_t * opt) const
 654 {
 655   //Print some relevant parameters set for the analysis
 656   if(! opt)
 657     return;
 658
 659   printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
 660   AliAnaPartCorrBaseClass::Print(" ");
 661
 662   printf("Select Calorimeter %s \n",fCalorimeter.Data());
 663   printf("Plots style macro %s \n",fStyleMacro.Data());
 664 }
 665
 666 //__________________________________________________________________
 667 void  AliAnaCalorimeterQA::MakeAnalysisFillHistograms()
 668 {
 669         //Do analysis and fill histograms
 670         TLorentzVector mom ;
 671         TLorentzVector mom2 ;
 672         //Play with the MC stack if available
 673         AliStack * stack = 0x0;
 674         TParticle * primary = 0x0;
 675         if(IsDataMC()) stack =  GetMCStack() ;
 676
 677         //Get List with clusters
 678         TRefArray * partList = new TRefArray();
 679         if(fCalorimeter == "EMCAL") partList = GetAODEMCAL();
 680         else if(fCalorimeter == "PHOS") partList = GetAODPHOS();
 681         else {
 682                 printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - Wrong calorimeter name <%s>, END\n", fCalorimeter.Data());
 683                 abort();
 684         }
 685
 686         if(!partList || partList->GetEntriesFast() == 0) return ;
 687
 688         Int_t nclusters = partList->GetEntriesFast() ;
 689         fhNClusters->Fill(nclusters);
 690
 691         if(GetDebug() > 0)
 692                 printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - In %s there are %d clusters \n", fCalorimeter.Data(), nclusters);
 693
 694         //Get vertex for photon momentum calculation
 695         Double_t v[3] ; //vertex ;
 696         GetReader()->GetVertex(v);
 697
 698         for(Int_t iclus = 0; iclus < partList->GetEntriesFast(); iclus++){
 699
 700                 AliAODCaloCluster * calo =  (AliAODCaloCluster*) (partList->At(iclus));
 701                 //if(calo->GetNCells() <= 2) continue;
 702                 //Get cluster kinematics
 703                 calo->GetMomentum(mom,v);
 704                 Float_t e   = mom.E();
 705                 //if(e < 0.5) continue;
 706                 //printf("e %2.2f, n %f\n",e, calo->GetNCells());
 707                 Float_t pt  = mom.Pt();
 708                 Float_t eta = mom.Eta();
 709                 Float_t phi = mom.Phi();
 710                 if(phi < 0) phi +=TMath::TwoPi();
 711
 712                 fhE     ->Fill(e);
 713                 fhPt    ->Fill(pt);
 714                 fhPhi   ->Fill(phi);
 715                 fhEta   ->Fill(eta);
 716                 fhEtaPhi->Fill(eta,phi);
 717
 718                 //matched cluster with tracks
 719                 if(calo->GetNTracksMatched() > 0){
 720                         fhECharged     ->Fill(e);
 721                         fhPtCharged    ->Fill(pt);
 722                         fhPhiCharged   ->Fill(phi);
 723                         fhEtaCharged   ->Fill(eta);
 724                         fhEtaPhiCharged->Fill(eta,phi);
 725                         if((!strcmp(GetReader()->GetInputEvent()->GetName(),"AliESDEvent"))) {
 726                                 AliESDEvent *esd = (AliESDEvent*) GetReader()->GetInputEvent();
 727                                 AliESDCaloCluster * esdcalo = (AliESDCaloCluster*) esd->GetCaloCluster(calo->GetID());
 728                                 Int_t trackIndex = esdcalo->GetTrackMatched();
 729                                 //printf("track index %d ntracks %d\n", trackIndex, esd->GetNumberOfTracks());
 730                                 if(trackIndex >= 0){
 731                                         AliESDtrack* track = (AliESDtrack*) esd->GetTrack(trackIndex);
 732                                         if (track && track->GetOuterParam() ) {
 733
 734                                                 Double_t tphi = track->GetOuterParam()->Phi();
 735                                                 Double_t teta = track->GetOuterParam()->Eta();
 736                                                 Double_t tmom = track->GetOuterParam()->P();
 737                                                 Double_t deta = teta - eta;
 738                                                 Double_t dphi = tphi - phi;
 739                                                 if(dphi > TMath::Pi()) dphi -= 2*TMath::Pi();
 740                                                 if(dphi < -TMath::Pi()) dphi += 2*TMath::Pi();
 741                                                 Double_t dR = sqrt(dphi*dphi + deta*deta);
 742
 743                                                 Double_t pOverE = tmom/e;
 744
 745                                                 fh1pOverE->Fill(pOverE);
 746                                                 fh1dR->Fill(dR);
 747                                                 fh2MatchdEdx->Fill(track->P(),track->GetTPCsignal());
 748
 749                                                 int nITS = track->GetNcls(0);
 750                                                 int nTPC = track->GetNcls(1);
 751                                                 if(dR < 0.02 && pOverE > 0.5 && pOverE < 1.5
 752                                                    && calo->GetNCells() > 1 && nITS > 3 && nTPC > 20) {
 753                                                         fh2EledEdx->Fill(track->P(),track->GetTPCsignal());
 754                                                         //cout<<track->P()<<" "<<track->P()<<endl;
 755                                                 }
 756                                         }
 757                                         else if(!track->GetOuterParam()){
 758                                           ULong_t status=AliESDtrack::kTPCrefit;
 759                                           status|=AliESDtrack::kITSrefit;
 760                                           //printf("track status %d\n", track->GetStatus() );
 761                                           fhEChargedNoOut     ->Fill(e);
 762                                           fhPtChargedNoOut     ->Fill(pt);
 763                                           fhPhiChargedNoOut    ->Fill(phi);
 764                                           fhEtaChargedNoOut    ->Fill(eta);
 765                                           fhEtaPhiChargedNoOut ->Fill(eta,phi);
 766                                           if ( (track->GetStatus() & status) == status) printf("ITS+TPC\n");
 767                                         }
 768                                         else {
 769                                                 Printf("ERROR: Could not receive track %d", trackIndex);
 770                                         }
 771                                 }// non negative track index
 772                         }//do only if input are ESDs
 773                 }// at least one track matched
 774
 775                 if(IsDataMC()){
 776                 //Play with the MC stack if available
 777                     Int_t label = calo->GetLabel(0);
 778
 779                 if(label < 0 || !stack) {
 780                                 printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() *** bad label or no stack ***:  label %d \n", label);
 781                                 continue;
 782                 }
 783
 784                 if(label >=  stack->GetNtrack()) {
 785                                 printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() *** large label ***:  label %d, n tracks %d \n", label, stack->GetNtrack());
 786                                 continue ;
 787                 }
 788                 //cout<<"LABEL > "<<label<<endl;
 789                 primary = GetMCStack()->Particle(label);
 790                     Int_t pdg = primary->GetPdgCode();
 791                     Float_t vx = primary->Vx();
 792                     Float_t vy = primary->Vy();
 793                     //Float_t vz = primary->Vz();
 794                     Float_t r = TMath::Sqrt(vx*vx + vy*vy);
 795                     if((pdg == 22 || TMath::Abs(pdg)==11) && primary->GetStatusCode()!=1) {
 796                       fhEMVxyz   ->Fill(vx,vy);//,vz);
 797                       fhEMR      ->Fill(e,r);
 798                     }
 799
 800                     //                  printf("reco e %f, pt %f, phi %f, eta %f \n", e, pt, phi, eta);
 801                     //                  printf("prim e %f, pt %f, phi %f, eta %f \n", primary->Energy(),primary->Pt() ,primary->Phi() ,primary->Eta() );
 802                     //                  printf("vertex: vx %f, vy %f, vz %f, r %f \n", vx, vy, vz, r);
 803
 804                     //Get final particle, no conversion products
 805                     Int_t status =  primary->GetStatusCode();
 806                     Int_t mother= primary->GetFirstMother();
 807                     Int_t finallabel = -1;
 808                     if(status == 0){
 809                       while(mother >= 0){
 810                         primary = GetMCStack()->Particle(mother);
 811                         finallabel = mother;
 812                         status = primary->GetStatusCode();
 813                         mother= primary->GetFirstMother();
 814                         if(status == 1) break;
 815                         //printf("mother %d\n",mother);
 816                       }
 817                     }
 818
 819                         fh2E      ->Fill(e, primary->Energy());
 820                         fh2Pt     ->Fill(pt, primary->Pt());
 821                         fh2Phi    ->Fill(phi, primary->Phi());
 822                         fh2Eta    ->Fill(eta, primary->Eta());
 823                         fhDeltaE  ->Fill(primary->Energy()-e);
 824                         fhDeltaPt ->Fill(primary->Pt()-pt);
 825                         fhDeltaPhi->Fill(primary->Phi()-phi);
 826                         fhDeltaEta->Fill(primary->Eta()-eta);
 827                         if(primary->Energy() > 0) fhRatioE  ->Fill(e/primary->Energy());
 828                     if(primary->Pt()     > 0) fhRatioPt ->Fill(pt/primary->Pt());
 829                     if(primary->Phi()    > 0) fhRatioPhi->Fill(phi/primary->Phi());
 830                     if(primary->Eta()    > 0) fhRatioEta->Fill(eta/primary->Eta());
 831
 832                     //cout<<"Final Label "<<finallabel<<" mother "<<mother<<endl;
 833                     pdg = primary->GetPdgCode();
 834                         Double_t  charge = TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
 835
 836                     if(pdg == 22){
 837                       //cout<<"pdg "<<pdg<<" status "<<status<<" "<<primary->GetStatusCode()<<endl;
 838                       TParticle *pi0 = 0x0;
 839                       TParticle *p1 = 0x0;
 840                       TParticle *p2 = 0x0;
 841                       TParticle * tmp = 0x0;
 842                       Int_t pdgpi0 = 0;
 843                       Int_t mothertmp = 0;
 844                       //Check if it is a decay photon from pi0, in this case, check if both decay contribute cluster
 845                       Bool_t ok1 = kFALSE;
 846                       Bool_t ok2 = kFALSE;
 847
 848                       if(calo->GetNLabel() > 1){
 849                         if(pdg !=111){
 850                           while(mother >= 0){
 851                             pi0 = GetMCStack()->Particle(mother);
 852                             pdgpi0 = pi0->GetPdgCode();
 853                             if(pdgpi0 == 111) break;
 854                             mother= pi0->GetFirstMother();
 855                             //printf("mother %d\n",mother);
 856                           }
 857                         }
 858                         else   pi0 = primary;
 859                         //cout<<"MOTHER PI0 LABEL "<<mother<<" pt" << pi0->Pt()<<endl;
 860                         if(pi0->GetNDaughters() == 2){
 861                           //  cout<<"pi0, 2 daughters "<<endl;
 862                           Int_t id1 = pi0->GetFirstDaughter();
 863                           Int_t id2 = pi0->GetFirstDaughter()+1;
 864                           p1=GetMCStack()->Particle(id1);
 865                           p2=GetMCStack()->Particle(id2);
 866
 867                           if(p1->GetFirstMother()!=p2->GetFirstMother()) cout <<"Decay photon mothers are not the same!!"<<endl;
 868                           if(p1->GetPdgCode()==22 && p2->GetPdgCode()==22){
 869                             // cout<<"2 photons, labels "<< id1<<" "<<id2<<endl;
 870                             for(UInt_t ilabel = 0; ilabel < calo->GetNLabel(); ilabel++){
 871                               Int_t iprim = calo->GetLabel(ilabel);
 872                               //cout<<"iprim "<<iprim<<endl;
 873                               if (iprim == id1) ok1 = kTRUE;
 874                               else if (iprim == id2) ok2 = kTRUE;
 875                               mothertmp = iprim;
 876                               while(mothertmp >= 0){
 877                                 tmp = GetMCStack()->Particle(mothertmp);
 878                                 mothertmp= tmp->GetFirstMother();
 879                                 //      cout<<"mothertmp "<<mothertmp<<" "<<tmp->GetName()<< " pt "<<tmp->Pt()<<endl;
 880                                 if (mothertmp == id1) ok1 = kTRUE;
 881                                 else if (mothertmp == id2) ok2 = kTRUE;
 882                                 if(ok1 && ok2) break;
 883                                 //printf("mother %d\n",mother);
 884                               }
 885                             }
 886                           }//2 photon daughters
 887                         }////mother pi0, 2 daughers
 888                       }//more than one contribution to clust
 889
 890                       if(ok1 && ok2){
 891                         //cout<<"Fill pi0"<< "E  "<< e <<" prim E "<<primary->Energy()<<endl;
 892                         fhPi0E     ->Fill(e,primary->Energy());
 893                         fhPi0Pt    ->Fill(pt,primary->Pt());
 894                         fhPi0Eta   ->Fill(eta,primary->Eta());
 895                         fhPi0Phi   ->Fill(phi,primary->Phi());
 896                       }
 897                       else{
 898                         fhGamE     ->Fill(e,primary->Energy());
 899                         fhGamPt    ->Fill(pt,primary->Pt());
 900                         fhGamEta   ->Fill(eta,primary->Eta());
 901                         fhGamPhi   ->Fill(phi,primary->Phi());
 902                         fhGamDeltaE  ->Fill(primary->Energy()-e);
 903                         fhGamDeltaPt ->Fill(primary->Pt()-pt);
 904                         fhGamDeltaPhi->Fill(primary->Phi()-phi);
 905                         fhGamDeltaEta->Fill(primary->Eta()-eta);
 906                         if(primary->Energy() > 0) fhGamRatioE  ->Fill(e/primary->Energy());
 907                         if(primary->Pt()     > 0) fhGamRatioPt ->Fill(pt/primary->Pt());
 908                         if(primary->Phi()    > 0) fhGamRatioPhi->Fill(phi/primary->Phi());
 909                         if(primary->Eta()    > 0) fhGamRatioEta->Fill(eta/primary->Eta());
 910                         fhGam2E      ->Fill(e, primary->Energy());
 911                         fhGam2Pt     ->Fill(pt, primary->Pt());
 912                         fhGam2Phi    ->Fill(phi, primary->Phi());
 913                         fhGam2Eta    ->Fill(eta, primary->Eta());
 914                       }
 915                     }//pdg == 22
 916                     else if(TMath::Abs(pdg) == 11) {
 917                       fhEleE     ->Fill(e,primary->Energy());
 918                       fhElePt    ->Fill(pt,primary->Pt());
 919                       fhEleEta   ->Fill(eta,primary->Eta());
 920                       fhElePhi   ->Fill(phi,primary->Phi());
 921                       fhEMVxyz   ->Fill(vx,vy);//,vz);
 922                       fhEMR      ->Fill(e,r);
 923                     }
 924                     else if(pdg == 111) {
 925                       fhPi0E     ->Fill(e,primary->Energy());
 926                       fhPi0Pt    ->Fill(pt,primary->Pt());
 927                       fhPi0Eta   ->Fill(eta,primary->Eta());
 928                       fhPi0Phi   ->Fill(phi,primary->Phi());
 929                     }
 930                     else if(charge == 0){
 931                       fhNeHadE     ->Fill(e,primary->Energy());
 932                       fhNeHadPt    ->Fill(pt,primary->Pt());
 933                       fhNeHadEta   ->Fill(eta,primary->Eta());
 934                       fhNeHadPhi   ->Fill(phi,primary->Phi());
 935                       fhHaVxyz     ->Fill(vx,vy);//,vz);
 936                       fhHaR        ->Fill(e,r);
 937                     }
 938                     else if(charge!=0){
 939                       fhChHadE     ->Fill(e,primary->Energy());
 940                       fhChHadPt    ->Fill(pt,primary->Pt());
 941                       fhChHadEta   ->Fill(eta,primary->Eta());
 942                       fhChHadPhi   ->Fill(phi,primary->Phi());
 943                       fhHaVxyz     ->Fill(vx,vy);//,vz);
 944                       fhHaR        ->Fill(e,r);
 945                     }
 946                 }//Work with stack also
 947
 948                 //Invariant mass
 949                 if (nclusters > 1 ) {
 950                         for(Int_t jclus = iclus + 1 ; jclus < nclusters ; jclus++) {
 951                                 AliAODCaloCluster * calo2 =  (AliAODCaloCluster*) (partList->At(jclus));
 952                                 //Get cluster kinematics
 953                                 calo2->GetMomentum(mom2,v);
 954                                 fhIM  ->Fill((mom+mom2).E(),(mom+mom2).M());
 955                                 fhAsym->Fill((mom+mom2).E(),TMath::Abs((e-mom2.E())/(e+mom2.E())));
 956                         }// 2nd cluster loop
 957                 }////more than 1 cluster in calorimeter
 958
 959                 //Cells per cluster
 960                 fhNCellsPerCluster->Fill(e,calo->GetNCells());
 961
 962         }// 1st cluster loop
 963
 964         //CaloCells
 965         AliAODCaloCells * cell = new AliAODCaloCells ;
 966         if(fCalorimeter == "PHOS")
 967                 cell = (AliAODCaloCells *) GetPHOSCells();
 968         else
 969                 cell = (AliAODCaloCells *) GetEMCALCells();
 970
 971         if(!cell) {
 972                 printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - ABORT: No CELLS available for analysis");
 973                 abort();
 974         }
 975
 976         //Some prints
 977         if(GetDebug() > 0 && cell )
 978                 printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - In ESD %s cell entries %d\n", fCalorimeter.Data(), cell->GetNumberOfCells());
 979
 980         Int_t ncells = cell->GetNumberOfCells() ;
 981         fhNCells->Fill(ncells) ;
 982
 983         for (Int_t iCell = 0; iCell < ncells; iCell++) {
 984                 if(GetDebug() > 2)  printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - Cell : amp %f, absId %d \n", cell->GetAmplitude(iCell), cell->GetCellNumber(iCell));
 985                 fhAmplitude->Fill(cell->GetAmplitude(iCell));
 986         }
 987
 988         if(GetDebug() > 0)
 989                 printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - end \n");
 990
 991         //Monte Carlo
 992         if(IsDataMC()){
 993                 //Play with the MC stack if available
 994                 for(Int_t i=8 ; i<stack->GetNprimary(); i++){
 995                         primary = stack->Particle(i) ;
 996
 997                         //if (!primary->IsPrimary()) continue;
 998                         if (TMath::Abs(primary->Eta()) > 1) continue;
 999
1000                         Int_t kf = primary->GetPdgCode();
1001                         //printf("kf %d\n",kf);
1002
1003                         Bool_t in = kTRUE;
1004                         primary->Momentum(mom);
1005                         if(IsFidutialCutOn()) in =  GetFidutialCut()->IsInFidutialCut(mom,fCalorimeter) ;
1006
1007                         if (kf==22) {
1008                                 fhGenGamPt ->Fill(primary->Pt());
1009                                 fhGenGamEta->Fill(primary->Eta());
1010                                 fhGenGamPhi->Fill(primary->Phi());
1011                                 if(in){
1012                                         fhGenGamAccE  ->Fill(primary->Energy());
1013                                         fhGenGamAccPt ->Fill(primary->Pt());
1014                                         fhGenGamAccEta->Fill(primary->Eta());
1015                                         fhGenGamAccPhi->Fill(primary->Phi());
1016                                 }
1017                         }
1018                         else if (kf==111) {
1019                                 fhGenPi0Pt ->Fill(primary->Pt());
1020                                 fhGenPi0Eta->Fill(primary->Eta());
1021                                 fhGenPi0Phi->Fill(primary->Phi());
1022                                 if(in){
1023                                         fhGenPi0AccE  ->Fill(primary->Energy());
1024                                         fhGenPi0AccPt ->Fill(primary->Pt());
1025                                         fhGenPi0AccEta->Fill(primary->Eta());
1026                                         fhGenPi0AccPhi->Fill(primary->Phi());
1027                                 }
1028                         }
1029                         else if (kf==221) {
1030                                 fhGenEtaPt ->Fill(primary->Pt());
1031                                 fhGenEtaEta->Fill(primary->Eta());
1032                                 fhGenEtaPhi->Fill(primary->Phi());
1033                         }
1034                         else if (kf==223) {
1035                                 fhGenOmegaPt ->Fill(primary->Pt());
1036                                 fhGenOmegaEta->Fill(primary->Eta());
1037                                 fhGenOmegaPhi->Fill(primary->Phi());
1038                         }
1039                         else if (TMath::Abs(kf)==11) {
1040                                 fhGenElePt ->Fill(primary->Pt());
1041                                 fhGenEleEta->Fill(primary->Eta());
1042                                 fhGenElePhi->Fill(primary->Phi());
1043                         }
1044                 } //primary loop
1045         } //Is data MC
1046
1047 }
1048
1049 //________________________________________________________________________
1050 void AliAnaCalorimeterQA::ReadHistograms(TList* outputList)
1051 {
1052         // Needed when Terminate is executed in distributed environment
1053         // Refill analysis histograms of this class with corresponding histograms in output list.
1054
1055         // Histograms of this analsys are kept in the same list as other analysis, recover the position of
1056         // the first one and then add the next
1057         Int_t index = outputList->IndexOf(outputList->FindObject(GetAddedHistogramsStringToName()+"hE"));
1058         //printf("Calo: %s, index: %d\n",fCalorimeter.Data(),index);
1059         //Read histograms, must be in the same order as in GetCreateOutputObject.
1060         fhE      = (TH1F *) outputList->At(index++);
1061         fhPt     = (TH1F *) outputList->At(index++);
1062         fhPhi    = (TH1F *) outputList->At(index++);
1063         fhEta    = (TH1F *) outputList->At(index++);
1064         fhEtaPhi = (TH2F *) outputList->At(index++);
1065
1066         fhECharged      = (TH1F *) outputList->At(index++);
1067         fhPtCharged     = (TH1F *) outputList->At(index++);
1068         fhPhiCharged    = (TH1F *) outputList->At(index++);
1069         fhEtaCharged    = (TH1F *) outputList->At(index++);
1070         fhEtaPhiCharged = (TH2F *) outputList->At(index++);
1071
1072         fhEChargedNoOut      = (TH1F *) outputList->At(index++);
1073         fhPtChargedNoOut     = (TH1F *) outputList->At(index++);
1074         fhPhiChargedNoOut    = (TH1F *) outputList->At(index++);
1075         fhEtaChargedNoOut    = (TH1F *) outputList->At(index++);
1076         fhEtaPhiChargedNoOut = (TH2F *) outputList->At(index++);
1077
1078         fhIM     = (TH2F *) outputList->At(index++);
1079         fhAsym   = (TH2F *) outputList->At(index++);
1080
1081         fhNCellsPerCluster = (TH2F *) outputList->At(index++);
1082         fhNClusters  = (TH1F *) outputList->At(index++);
1083         fhNCells     = (TH1F *) outputList->At(index++);
1084         fhAmplitude  = (TH1F *) outputList->At(index++);
1085
1086         if(IsDataMC()){
1087                 fhDeltaE   = (TH1F *) outputList->At(index++);
1088                 fhDeltaPt  = (TH1F *) outputList->At(index++);
1089                 fhDeltaPhi = (TH1F *) outputList->At(index++);
1090                 fhDeltaEta = (TH1F *) outputList->At(index++);
1091
1092                 fhRatioE   = (TH1F *) outputList->At(index++);
1093                 fhRatioPt  = (TH1F *) outputList->At(index++);
1094                 fhRatioPhi = (TH1F *) outputList->At(index++);
1095                 fhRatioEta = (TH1F *) outputList->At(index++);
1096
1097                 fh2E       = (TH2F *) outputList->At(index++);
1098                 fh2Pt      = (TH2F *) outputList->At(index++);
1099                 fh2Phi     = (TH2F *) outputList->At(index++);
1100                 fh2Eta     = (TH2F *) outputList->At(index++);
1101
1102                 fhGamE     = (TH2F *) outputList->At(index++);
1103                 fhGamPt    = (TH2F *) outputList->At(index++);
1104                 fhGamPhi   = (TH2F *) outputList->At(index++);
1105                 fhGamEta   = (TH2F *) outputList->At(index++);
1106
1107                 fhGamDeltaE   = (TH1F *) outputList->At(index++);
1108                 fhGamDeltaPt  = (TH1F *) outputList->At(index++);
1109                 fhGamDeltaPhi = (TH1F *) outputList->At(index++);
1110                 fhGamDeltaEta = (TH1F *) outputList->At(index++);
1111
1112                 fhGamRatioE   = (TH1F *) outputList->At(index++);
1113                 fhGamRatioPt  = (TH1F *) outputList->At(index++);
1114                 fhGamRatioPhi = (TH1F *) outputList->At(index++);
1115                 fhGamRatioEta = (TH1F *) outputList->At(index++);
1116
1117                 fhGam2E       = (TH2F *) outputList->At(index++);
1118                 fhGam2Pt      = (TH2F *) outputList->At(index++);
1119                 fhGam2Phi     = (TH2F *) outputList->At(index++);
1120                 fhGam2Eta     = (TH2F *) outputList->At(index++);
1121
1122                 fhPi0E     = (TH2F *) outputList->At(index++);
1123                 fhPi0Pt    = (TH2F *) outputList->At(index++);
1124                 fhPi0Phi   = (TH2F *) outputList->At(index++);
1125                 fhPi0Eta   = (TH2F *) outputList->At(index++);
1126
1127                 fhEleE     = (TH2F *) outputList->At(index++);
1128                 fhElePt    = (TH2F *) outputList->At(index++);
1129                 fhElePhi   = (TH2F *) outputList->At(index++);
1130                 fhEleEta   = (TH2F *) outputList->At(index++);
1131
1132                 fhNeHadE     = (TH2F *) outputList->At(index++);
1133                 fhNeHadPt    = (TH2F *) outputList->At(index++);
1134                 fhNeHadPhi   = (TH2F *) outputList->At(index++);
1135                 fhNeHadEta   = (TH2F *) outputList->At(index++);
1136
1137                 fhChHadE     = (TH2F *) outputList->At(index++);
1138                 fhChHadPt    = (TH2F *) outputList->At(index++);
1139                 fhChHadPhi   = (TH2F *) outputList->At(index++);
1140                 fhChHadEta   = (TH2F *) outputList->At(index++);
1141
1142 //              fhEMVxyz     = (TH3F *) outputList->At(index++);
1143 //              fhHaVxyz     = (TH3F *) outputList->At(index++);
1144
1145                 fhEMVxyz     = (TH2F *) outputList->At(index++);
1146                 fhHaVxyz     = (TH2F *) outputList->At(index++);
1147                 fhEMR        = (TH2F *) outputList->At(index++);
1148                 fhHaR        = (TH2F *) outputList->At(index++);
1149
1150                 fhGenGamPt    = (TH1F *) outputList->At(index++);
1151                 fhGenGamEta   = (TH1F *) outputList->At(index++);
1152                 fhGenGamPhi   = (TH1F *) outputList->At(index++);
1153
1154                 fhGenPi0Pt    = (TH1F *) outputList->At(index++);
1155                 fhGenPi0Eta   = (TH1F *) outputList->At(index++);
1156                 fhGenPi0Phi   = (TH1F *) outputList->At(index++);
1157
1158                 fhGenEtaPt    = (TH1F *) outputList->At(index++);
1159                 fhGenEtaEta   = (TH1F *) outputList->At(index++);
1160                 fhGenEtaPhi   = (TH1F *) outputList->At(index++);
1161
1162                 fhGenOmegaPt  = (TH1F *) outputList->At(index++);
1163                 fhGenOmegaEta = (TH1F *) outputList->At(index++);
1164                 fhGenOmegaPhi = (TH1F *) outputList->At(index++);
1165
1166                 fhGenElePt    = (TH1F *) outputList->At(index++);
1167                 fhGenEleEta   = (TH1F *) outputList->At(index++);
1168                 fhGenElePhi   = (TH1F *) outputList->At(index++);
1169
1170                 fhGenGamAccE   = (TH1F *) outputList->At(index++);
1171                 fhGenGamAccPt  = (TH1F *) outputList->At(index++);
1172                 fhGenGamAccEta = (TH1F *) outputList->At(index++);
1173                 fhGenGamAccPhi = (TH1F *) outputList->At(index++);
1174
1175                 fhGenPi0AccE   = (TH1F *) outputList->At(index++);
1176                 fhGenPi0AccPt  = (TH1F *) outputList->At(index++);
1177                 fhGenPi0AccEta = (TH1F *) outputList->At(index++);
1178                 fhGenPi0AccPhi = (TH1F *) outputList->At(index++);
1179
1180         }//Is data MC
1181
1182         fh1pOverE =    (TH1F *) outputList->At(index++);
1183         fh1dR =        (TH1F *) outputList->At(index++);
1184         fh2MatchdEdx = (TH2F *) outputList->At(index++);
1185         fh2EledEdx =   (TH2F *) outputList->At(index++);
1186
1187 //      for(Int_t i = 0;  i<index ; i++) cout<<outputList->At(i)->GetName()<<endl;
1188 }
1189
1190 //__________________________________________________________________
1191 void  AliAnaCalorimeterQA::Terminate(TList* outputList)
1192 {
1193
1194         //Do some plots to end
1195          if(fStyleMacro!="")gROOT->Macro(fStyleMacro);
1196         //Recover histograms from output histograms list, needed for distributed analysis.
1197         ReadHistograms(outputList);
1198
1199         //printf(" AliAnaCalorimeterQA::Terminate()  *** %s Report:", GetName()) ;
1200         //printf(" AliAnaCalorimeterQA::Terminate()        pt         : %5.3f , RMS : %5.3f \n", fhPt->GetMean(),   fhPt->GetRMS() ) ;
1201
1202         char name[128];
1203         char cname[128];
1204
1205         //Reconstructed distributions
1206         //printf("c\n");
1207         sprintf(cname,"QA_%s_rec",fCalorimeter.Data());
1208         TCanvas  * c = new TCanvas(cname, "Reconstructed distributions", 400, 400) ;
1209         c->Divide(2, 2);
1210
1211         c->cd(1) ;
1212         gPad->SetLogy();
1213         fhE->SetLineColor(4);
1214         fhE->Draw();
1215
1216         c->cd(2) ;
1217         gPad->SetLogy();
1218         fhPt->SetLineColor(4);
1219         fhPt->Draw();
1220
1221         c->cd(3) ;
1222         fhPhi->SetLineColor(4);
1223         fhPhi->Draw();
1224
1225         c->cd(4) ;
1226         fhEta->SetLineColor(4);
1227         fhEta->Draw();
1228
1229         sprintf(name,"QA_%s_ReconstructedDistributions.eps",fCalorimeter.Data());
1230         c->Print(name);
1231
1232         //Reconstructed distributions, matched with tracks
1233         //printf("c2\n");
1234         sprintf(cname,"QA_%s_rectrackmatch",fCalorimeter.Data());
1235         TCanvas  * c2 = new TCanvas(cname, "Reconstructed distributions, matched with tracks", 400, 400) ;
1236         c2->Divide(2, 2);
1237
1238         c2->cd(1) ;
1239         gPad->SetLogy();
1240         fhECharged->SetLineColor(4);
1241         fhECharged->Draw();
1242
1243         c2->cd(2) ;
1244         gPad->SetLogy();
1245         fhPtCharged->SetLineColor(4);
1246         fhPtCharged->Draw();
1247
1248         c2->cd(3) ;
1249         fhPhiCharged->SetLineColor(4);
1250         fhPhiCharged->Draw();
1251
1252         c2->cd(4) ;
1253         fhEtaCharged->SetLineColor(4);
1254         fhEtaCharged->Draw();
1255
1256         sprintf(name,"QA_%s_ReconstructedDistributions_TrackMatched.eps",fCalorimeter.Data());
1257         c2->Print(name);
1258
1259         TH1F *  hEChargedClone   = (TH1F*)   fhECharged->Clone("EChargedClone");
1260         TH1F *  hPtChargedClone  = (TH1F*)   fhPtCharged->Clone("PtChargedClone");
1261         TH1F *  hEtaChargedClone = (TH1F*)   fhEtaCharged->Clone("EtaChargedClone");
1262         TH1F *  hPhiChargedClone = (TH1F*)   fhPhiCharged->Clone("PhiChargedClone");
1263
1264         TH1F *  hEChargedClone2   = (TH1F*)   fhECharged->Clone("EChargedClone2");
1265         TH1F *  hPtChargedClone2  = (TH1F*)   fhPtCharged->Clone("PtChargedClone2");
1266         TH1F *  hEtaChargedClone2 = (TH1F*)   fhEtaCharged->Clone("EtaChargedClone2");
1267         TH1F *  hPhiChargedClone2 = (TH1F*)   fhPhiCharged->Clone("PhiChargedClone2");
1268
1269         //Ratio: reconstructed track matched/ all reconstructed
1270         //printf("c3\n");
1271         sprintf(cname,"QA_%s_rectrackmatchrat",fCalorimeter.Data());
1272         TCanvas  * c3 = new TCanvas(cname, "Ratio: reconstructed track matched/ all reconstructed", 400, 400) ;
1273         c3->Divide(2, 2);
1274
1275         c3->cd(1) ;
1276         hEChargedClone->SetYTitle("track matched / all");
1277         hEChargedClone->Divide(fhE);
1278         hEChargedClone->Draw();
1279
1280         c3->cd(2) ;
1281         hPtChargedClone->SetYTitle("track matched / all");
1282         hPtChargedClone->Divide(fhPt);
1283         hPtChargedClone->Draw();
1284
1285         c3->cd(3) ;
1286         hPhiChargedClone->SetYTitle("track matched / all");
1287         hPhiChargedClone->Divide(fhPhi);
1288         hPhiChargedClone->Draw();
1289
1290         c3->cd(4) ;
1291         hEtaChargedClone->SetYTitle("track matched / all");
1292         hEtaChargedClone->Divide(fhEta);
1293         hEtaChargedClone->Draw();
1294
1295         sprintf(name,"QA_%s_RatioReconstructedMatchedDistributions.eps",fCalorimeter.Data());
1296         c3->Print(name);
1297
1298         //Ratio: reconstructed track matched (minus no track param) / all
1299         //printf("c333\n");
1300         sprintf(cname,"QA_%s_rectrackmatchratout",fCalorimeter.Data());
1301         TCanvas  * c333 = new TCanvas(cname, "Ratio: reconstructed track matched (with outer track param)/ all", 400, 400) ;
1302         c333->Divide(2, 2);
1303
1304         c333->cd(1) ;
1305         hEChargedClone2->Add(fhEChargedNoOut,-1);
1306         hEChargedClone2->SetYTitle("track matched / all");
1307         hEChargedClone2->Divide(fhE);
1308         hEChargedClone2->Draw();
1309
1310         c333->cd(2) ;
1311         hPtChargedClone2->Add(fhPtChargedNoOut,-1);
1312         hPtChargedClone2->SetYTitle("track matched / all");
1313         hPtChargedClone2->Divide(fhPt);
1314         hPtChargedClone2->Draw();
1315
1316         c333->cd(3) ;
1317         hPhiChargedClone2->Add(fhPhiChargedNoOut,-1);
1318         hPhiChargedClone2->SetYTitle("track matched / all");
1319         hPhiChargedClone2->Divide(fhPhi);
1320         hPhiChargedClone2->Draw();
1321
1322         c333->cd(4) ;
1323         hEtaChargedClone2->Add(fhEtaChargedNoOut,-1);
1324         hEtaChargedClone2->SetYTitle("track matched / all");
1325         hEtaChargedClone2->Divide(fhEta);
1326         hEtaChargedClone2->Draw();
1327
1328         sprintf(name,"QA_%s_RatioReconstructedMatchedDistributionsOuter.eps",fCalorimeter.Data());
1329         c333->Print(name);
1330
1331         //Reconstructed distributions, matched with tracks
1332         //printf("c2\n");
1333         sprintf(cname,"QA_%s_rectrackmatch_noout",fCalorimeter.Data());
1334         TCanvas  * c22 = new TCanvas(cname, "Reconstructed distributions, matched with tracks, no outer track param", 400, 400) ;
1335         c22->Divide(2, 2);
1336
1337         c22->cd(1) ;
1338         gPad->SetLogy();
1339         fhEChargedNoOut->SetLineColor(4);
1340         fhEChargedNoOut->Draw();
1341
1342         c22->cd(2) ;
1343         gPad->SetLogy();
1344         fhPtChargedNoOut->SetLineColor(4);
1345         fhPtChargedNoOut->Draw();
1346
1347         c22->cd(3) ;
1348         fhPhiChargedNoOut->SetLineColor(4);
1349         fhPhiChargedNoOut->Draw();
1350
1351         c22->cd(4) ;
1352         fhEtaChargedNoOut->SetLineColor(4);
1353         fhEtaChargedNoOut->Draw();
1354
1355         sprintf(name,"QA_%s_ReconstructedDistributions_TrackMatched_NoOutParam.eps",fCalorimeter.Data());
1356         c22->Print(name);
1357
1358         //Ratio: reconstructed track matched/ all reconstructed
1359         //printf("c33\n");
1360         TH1F *  hEChargedNoOutClone   = (TH1F*)   fhEChargedNoOut->Clone("EChargedNoOutClone");
1361         TH1F *  hPtChargedNoOutClone  = (TH1F*)   fhPtChargedNoOut->Clone("PtChargedNoOutClone");
1362         TH1F *  hEtaChargedNoOutClone = (TH1F*)   fhEtaChargedNoOut->Clone("EtaChargedNoOutClone");
1363         TH1F *  hPhiChargedNoOutClone = (TH1F*)   fhPhiChargedNoOut->Clone("PhiChargedNoOutClone");
1364
1365         sprintf(cname,"QA_%s_rectrackmatchratnoout",fCalorimeter.Data());
1366         TCanvas  * c33 = new TCanvas(cname, "Ratio: reconstructed track matched/ all reconstructed", 400, 400) ;
1367         c33->Divide(2, 2);
1368
1369         c33->cd(1) ;
1370         hEChargedNoOutClone->SetYTitle("track matched no out/ all matched");
1371         hEChargedNoOutClone->Divide(fhECharged);
1372         hEChargedNoOutClone->Draw();
1373
1374         c33->cd(2) ;
1375         hPtChargedNoOutClone->SetYTitle("track matched no out / all matched");
1376         hPtChargedNoOutClone->Divide(fhPtCharged);
1377         hPtChargedNoOutClone->Draw();
1378
1379         c33->cd(3) ;
1380         hPhiChargedNoOutClone->SetYTitle("track matched no out/ all matched");
1381         hPhiChargedNoOutClone->Divide(fhPhiCharged);
1382         hPhiChargedNoOutClone->Draw();
1383
1384         c33->cd(4) ;
1385         hEtaChargedNoOutClone->SetYTitle("track matched no out/ all matched");
1386         hEtaChargedNoOutClone->Divide(fhEtaCharged);
1387         hEtaChargedNoOutClone->Draw();
1388
1389         sprintf(name,"QA_%s_RatioMatchedDistributionsAllToNoOut.eps",fCalorimeter.Data());
1390         c33->Print(name);
1391
1392         //eta vs phi
1393         //printf("c4\n");
1394         sprintf(cname,"QA_%s_etavsphi",fCalorimeter.Data());
1395         TCanvas  * c4 = new TCanvas(cname, "reconstructed #eta vs #phi", 200, 600) ;
1396         c4->Divide(1, 3);
1397
1398         c4->cd(1) ;
1399         fhEtaPhi->Draw("cont");
1400
1401         c4->cd(2) ;
1402         fhEtaPhiCharged->Draw("cont");
1403
1404         c4->cd(3) ;
1405         fhEtaPhiChargedNoOut->Draw("cont");
1406
1407         sprintf(name,"QA_%s_ReconstructedEtaVsPhi.eps",fCalorimeter.Data());
1408         c4->Print(name);
1409
1410         //Invariant mass
1411         Int_t binmin = -1;
1412         Int_t binmax = -1;
1413
1414         if(fhIM->GetEntries() > 1){
1415                 Int_t nebins  = fhIM->GetNbinsX();
1416                 Int_t emax = (Int_t) fhIM->GetXaxis()->GetXmax();
1417                 Int_t emin = (Int_t) fhIM->GetXaxis()->GetXmin();
1418                 if (emin != 0 ) printf("emin != 0 \n");
1419                 //printf("IM: nBinsX %d, emin %2.2f, emax %2.2f\n",nebins,emin,emax);
1420
1421                 sprintf(cname,"QA_%s_IM",fCalorimeter.Data());
1422         //      printf("c5\n");
1423                 TCanvas  * c5 = new TCanvas(cname, "Invariant mass", 400, 400) ;
1424                 c5->Divide(2, 2);
1425
1426                 c5->cd(1) ;
1427                 fhIM->SetLineColor(4);
1428                 fhIM->Draw();
1429
1430                 c5->cd(2) ;
1431                 binmin = 0;
1432                 binmax =  (Int_t) (5-emin)*nebins/emax;
1433                 TH1D *pyim5 = fhIM->ProjectionY("pyim5",binmin,binmax);
1434                 pyim5->SetTitle("E_{pair} < 5 GeV");
1435                 pyim5->SetLineColor(4);
1436                 pyim5->Draw();
1437
1438                 c5->cd(3) ;
1439                 binmin =  (Int_t) (5-emin)*nebins/emax;
1440                 binmax =  (Int_t) (10-emin)*nebins/emax;
1441                 TH1D *pyim510 = fhIM->ProjectionY("pyim5_10",binmin,binmax);
1442                 pyim510->SetTitle("5 < E_{pair} < 10 GeV");
1443                 pyim510->SetLineColor(4);
1444                 pyim510->Draw();
1445
1446                 c5->cd(4) ;
1447                 binmin =  (Int_t) (10-emin)*nebins/emax;
1448                 binmax = -1;
1449                 TH1D *pyim10 = fhIM->ProjectionY("pyim10",binmin,binmax);
1450                 pyim10->SetTitle("E_{pair} > 10 GeV");
1451                 pyim10->SetLineColor(4);
1452                 pyim10->Draw();
1453
1454                 sprintf(name,"QA_%s_InvariantMass.eps",fCalorimeter.Data());
1455                 c5->Print(name);
1456         }
1457
1458         //Asymmetry
1459         if(fhAsym->GetEntries() > 1){
1460                 Int_t nebins  = fhAsym->GetNbinsX();
1461                 Int_t emax = (Int_t) fhAsym->GetXaxis()->GetXmax();
1462                 Int_t emin = (Int_t) fhAsym->GetXaxis()->GetXmin();
1463                 if (emin != 0 ) printf("emin != 0 \n");
1464                 //printf("Asym: nBinsX %d, emin %2.2f, emax %2.2f\n",nebins,emin,emax);
1465
1466                 sprintf(cname,"QA_%s_Asym",fCalorimeter.Data());
1467                 //      printf("c5\n");
1468                 TCanvas  * c5b = new TCanvas(cname, "Asymmetry", 400, 400) ;
1469                 c5b->Divide(2, 2);
1470
1471                 c5b->cd(1) ;
1472                 fhAsym->SetLineColor(4);
1473                 fhAsym->Draw();
1474
1475                 c5b->cd(2) ;
1476                 binmin = 0;
1477                 binmax = (Int_t) (5-emin)*nebins/emax;
1478                 TH1D *pyAsym5 = fhAsym->ProjectionY("pyAsym5",binmin,binmax);
1479                 pyAsym5->SetTitle("E_{pair} < 5 GeV");
1480                 pyAsym5->SetLineColor(4);
1481                 pyAsym5->Draw();
1482
1483                 c5b->cd(3) ;
1484                 binmin = (Int_t) (5-emin)*nebins/emax;
1485                 binmax = (Int_t) (10-emin)*nebins/emax;
1486                 TH1D *pyAsym510 = fhAsym->ProjectionY("pyAsym5_10",binmin,binmax);
1487                 pyAsym510->SetTitle("5 < E_{pair} < 10 GeV");
1488                 pyAsym510->SetLineColor(4);
1489                 pyAsym510->Draw();
1490
1491                 c5b->cd(4) ;
1492                 binmin = (Int_t) (10-emin)*nebins/emax;
1493                 binmax = -1;
1494                 TH1D *pyAsym10 = fhAsym->ProjectionY("pyAsym10",binmin,binmax);
1495                 pyAsym10->SetTitle("E_{pair} > 10 GeV");
1496                 pyAsym10->SetLineColor(4);
1497                 pyAsym10->Draw();
1498
1499                 sprintf(name,"QA_%s_Asymmetry.eps",fCalorimeter.Data());
1500                 c5b->Print(name);
1501         }
1502
1503
1504         //Reconstructed vs MC distributions
1505         //printf("c6\n");
1506         sprintf(cname,"QA_%s_recvsmc",fCalorimeter.Data());
1507         TCanvas  * c6 = new TCanvas(cname, "Reconstructed vs MC distributions", 400, 400) ;
1508         c6->Divide(2, 2);
1509
1510         c6->cd(1) ;
1511         fh2E->SetLineColor(4);
1512         fh2E->Draw();
1513
1514         c6->cd(2) ;
1515         fh2Pt->SetLineColor(4);
1516         fh2Pt->Draw();
1517
1518         c6->cd(3) ;
1519         fh2Phi->SetLineColor(4);
1520         fh2Phi->Draw();
1521
1522         c6->cd(4) ;
1523         fh2Eta->SetLineColor(4);
1524         fh2Eta->Draw();
1525
1526         sprintf(name,"QA_%s_ReconstructedVSMCDistributions.eps",fCalorimeter.Data());
1527         c6->Print(name);
1528
1529
1530         //Reconstructed vs MC distributions
1531         //printf("c6\n");
1532         sprintf(cname,"QA_%s_gamrecvsmc",fCalorimeter.Data());
1533         TCanvas  * c6Gam = new TCanvas(cname, "Reconstructed vs MC distributions", 400, 400) ;
1534         c6Gam->Divide(2, 2);
1535
1536         c6Gam->cd(1) ;
1537         fhGam2E->SetLineColor(4);
1538         fhGam2E->Draw();
1539
1540         c6Gam->cd(2) ;
1541         fhGam2Pt->SetLineColor(4);
1542         fhGam2Pt->Draw();
1543
1544         c6Gam->cd(3) ;
1545         fhGam2Phi->SetLineColor(4);
1546         fhGam2Phi->Draw();
1547
1548         c6Gam->cd(4) ;
1549         fhGam2Eta->SetLineColor(4);
1550         fhGam2Eta->Draw();
1551
1552         sprintf(name,"QA_%s_GammaReconstructedVSMCDistributions.eps",fCalorimeter.Data());
1553         c6->Print(name);
1554
1555         //Generated - reconstructed
1556         //printf("c7\n");
1557         sprintf(cname,"QA_%s_diffgenrec",fCalorimeter.Data());
1558         TCanvas  * c7 = new TCanvas(cname, "generated - reconstructed", 400, 400) ;
1559         c7->Divide(2, 2);
1560
1561         c7->cd(1) ;
1562         gPad->SetLogy();
1563         fhGamDeltaE->SetLineColor(4);
1564         fhDeltaE->Draw();
1565         fhGamDeltaE->Draw("same");
1566
1567         TLegend pLegendd(0.65,0.55,0.9,0.8);
1568         pLegendd.SetTextSize(0.06);
1569         pLegendd.AddEntry(fhDeltaE,"all","L");
1570         pLegendd.AddEntry(fhGamDeltaE,"from  #gamma","L");
1571         pLegendd.SetFillColor(10);
1572         pLegendd.SetBorderSize(1);
1573         pLegendd.Draw();
1574
1575         c7->cd(2) ;
1576         gPad->SetLogy();
1577         fhGamDeltaPt->SetLineColor(4);
1578         fhDeltaPt->Draw();
1579         fhGamDeltaPt->Draw("same");
1580
1581         c7->cd(3) ;
1582         fhGamDeltaPhi->SetLineColor(4);
1583         fhDeltaPhi->Draw();
1584         fhGamDeltaPhi->Draw("same");
1585
1586         c7->cd(4) ;
1587         fhGamDeltaEta->SetLineColor(4);
1588         fhDeltaEta->Draw();
1589         fhGamDeltaEta->Draw("same");
1590
1591         sprintf(name,"QA_%s_DiffGeneratedReconstructed.eps",fCalorimeter.Data());
1592         c7->Print(name);
1593
1594
1595         // Reconstructed / Generated
1596         //printf("c8\n");
1597         sprintf(cname,"QA_%s_ratiorecgen",fCalorimeter.Data());
1598         TCanvas  * c8 = new TCanvas(cname, " reconstructed / generated", 400, 400) ;
1599         c8->Divide(2, 2);
1600
1601         c8->cd(1) ;
1602         gPad->SetLogy();
1603         fhGamRatioE->SetLineColor(4);
1604         fhRatioE->Draw();
1605         fhGamRatioE->Draw("same");
1606
1607         TLegend pLegendr(0.65,0.55,0.9,0.8);
1608         pLegendr.SetTextSize(0.06);
1609         pLegendr.AddEntry(fhRatioE,"all","L");
1610         pLegendr.AddEntry(fhGamRatioE,"from  #gamma","L");
1611         pLegendr.SetFillColor(10);
1612         pLegendr.SetBorderSize(1);
1613         pLegendr.Draw();
1614
1615         c8->cd(2) ;
1616         gPad->SetLogy();
1617         fhGamRatioPt->SetLineColor(4);
1618         fhRatioPt->Draw();
1619         fhGamRatioPt->Draw("same");
1620
1621         c8->cd(3) ;
1622         fhGamRatioPhi->SetLineColor(4);
1623         fhRatioPhi->Draw();
1624         fhGamRatioPhi->Draw("same");
1625
1626         c8->cd(4) ;
1627         fhGamRatioEta->SetLineColor(4);
1628         fhRatioEta->Draw();
1629         fhGamRatioEta->Draw("same");
1630
1631         sprintf(name,"QA_%s_ReconstructedDivGenerated.eps",fCalorimeter.Data());
1632         c8->Print(name);
1633
1634         // CaloClusters CaloCells
1635         //printf("c9\n");
1636         sprintf(cname,"QA_%s_nclustercellsamp",fCalorimeter.Data());
1637         TCanvas  * c9 = new TCanvas(cname, " CaloClusters and CaloCells", 400, 400) ;
1638         c9->Divide(2, 2);
1639
1640         c9->cd(1) ;
1641         gPad->SetLogy();
1642         gPad->SetLogx();
1643         fhNClusters->SetLineColor(4);
1644         fhNClusters->Draw();
1645
1646         c9->cd(2) ;
1647         gPad->SetLogy();
1648         gPad->SetLogx();
1649         fhNCells->SetLineColor(4);
1650         fhNCells->Draw();
1651
1652         c9->cd(3) ;
1653         gPad->SetLogy();
1654         gPad->SetLogx();
1655         fhNCellsPerCluster->SetLineColor(4);
1656         fhNCellsPerCluster->Draw();
1657
1658         c9->cd(4) ;
1659         gPad->SetLogy();
1660         //gPad->SetLogx();
1661         fhAmplitude->SetLineColor(4);
1662         fhAmplitude->Draw();
1663
1664         sprintf(name,"QA_%s_CaloClustersAndCaloCells.eps",fCalorimeter.Data());
1665         c9->Print(name);
1666
1667
1668         //MC
1669
1670         //Generated distributions
1671         //printf("c1\n");
1672         sprintf(cname,"QA_%s_gen",fCalorimeter.Data());
1673         TCanvas  * c10 = new TCanvas(cname, "Generated distributions", 600, 200) ;
1674         c10->Divide(3, 1);
1675
1676         c10->cd(1) ;
1677         gPad->SetLogy();
1678         TH1F * haxispt  = (TH1F*) fhGenPi0Pt->Clone("axispt");
1679         haxispt->SetTitle("Generated Particles p_{T}, |#eta| < 1");
1680         fhGenPi0Pt->SetLineColor(1);
1681         fhGenGamPt->SetLineColor(4);
1682         fhGenEtaPt->SetLineColor(2);
1683         fhGenOmegaPt->SetLineColor(7);
1684         fhGenElePt->SetLineColor(6);
1685
1686         //Select the maximum of the histogram to show all lines.
1687         if(fhGenPi0Pt->GetMaximum() >= fhGenGamPt->GetMaximum() && fhGenPi0Pt->GetMaximum() >= fhGenEtaPt->GetMaximum() &&
1688            fhGenPi0Pt->GetMaximum() >= fhGenOmegaPt->GetMaximum() && fhGenPi0Pt->GetMaximum() >= fhGenElePt->GetMaximum())
1689                 haxispt->SetMaximum(fhGenPi0Pt->GetMaximum());
1690         else if(fhGenGamPt->GetMaximum() >= fhGenPi0Pt->GetMaximum() && fhGenGamPt->GetMaximum() >= fhGenEtaPt->GetMaximum() &&
1691            fhGenGamPt->GetMaximum() >= fhGenOmegaPt->GetMaximum() && fhGenGamPt->GetMaximum() >= fhGenElePt->GetMaximum())
1692                 haxispt->SetMaximum(fhGenGamPt->GetMaximum());
1693         else if(fhGenEtaPt->GetMaximum() >= fhGenPi0Pt->GetMaximum() && fhGenEtaPt->GetMaximum() >= fhGenGamPt->GetMaximum() &&
1694                         fhGenEtaPt->GetMaximum() >= fhGenOmegaPt->GetMaximum() && fhGenEtaPt->GetMaximum() >= fhGenElePt->GetMaximum())
1695                 haxispt->SetMaximum(fhGenEtaPt->GetMaximum());
1696         else if(fhGenOmegaPt->GetMaximum() >= fhGenPi0Pt->GetMaximum() && fhGenOmegaPt->GetMaximum() >= fhGenEtaPt->GetMaximum() &&
1697                         fhGenOmegaPt->GetMaximum() >= fhGenGamPt->GetMaximum() && fhGenOmegaPt->GetMaximum() >= fhGenElePt->GetMaximum())
1698                 haxispt->SetMaximum(fhGenOmegaPt->GetMaximum());
1699         else if(fhGenElePt->GetMaximum() >= fhGenPi0Pt->GetMaximum() && fhGenElePt->GetMaximum() >= fhGenEtaPt->GetMaximum() &&
1700                         fhGenElePt->GetMaximum() >= fhGenOmegaPt->GetMaximum() && fhGenElePt->GetMaximum() >= fhGenGamPt->GetMaximum())
1701                 haxispt->SetMaximum(fhGenElePt->GetMaximum());
1702
1703         haxispt->Draw("axis");
1704         fhGenPi0Pt->Draw("same");
1705         fhGenGamPt->Draw("same");
1706         fhGenEtaPt->Draw("same");
1707         fhGenOmegaPt->Draw("same");
1708         fhGenElePt->Draw("same");
1709
1710         TLegend pLegend(0.75,0.45,0.9,0.8);
1711         pLegend.SetTextSize(0.06);
1712         pLegend.AddEntry(fhGenPi0Pt,"#pi^{0}","L");
1713         pLegend.AddEntry(fhGenGamPt,"#gamma","L");
1714         pLegend.AddEntry(fhGenEtaPt,"#eta","L");
1715         pLegend.AddEntry(fhGenOmegaPt,"#omega","L");
1716         pLegend.AddEntry(fhGenElePt,"e^{#pm}","L");
1717         pLegend.SetFillColor(10);
1718         pLegend.SetBorderSize(1);
1719         pLegend.Draw();
1720
1721         c10->cd(2) ;
1722         gPad->SetLogy();
1723         TH1F * haxiseta  = (TH1F*) fhGenPi0Eta->Clone("axiseta");
1724         haxiseta->SetTitle("Generated Particles #eta, |#eta| < 1");
1725         fhGenPi0Eta->SetLineColor(1);
1726         fhGenGamEta->SetLineColor(4);
1727         fhGenEtaEta->SetLineColor(2);
1728         fhGenOmegaEta->SetLineColor(7);
1729         fhGenEleEta->SetLineColor(6);
1730         //Select the maximum of the histogram to show all lines.
1731         if(fhGenPi0Eta->GetMaximum() >= fhGenGamEta->GetMaximum() && fhGenPi0Eta->GetMaximum() >= fhGenEtaEta->GetMaximum() &&
1732            fhGenPi0Eta->GetMaximum() >= fhGenOmegaEta->GetMaximum() && fhGenPi0Eta->GetMaximum() >= fhGenEleEta->GetMaximum())
1733                 haxiseta->SetMaximum(fhGenPi0Eta->GetMaximum());
1734         else if(fhGenGamEta->GetMaximum() >= fhGenPi0Eta->GetMaximum() && fhGenGamEta->GetMaximum() >= fhGenEtaEta->GetMaximum() &&
1735                         fhGenGamEta->GetMaximum() >= fhGenOmegaEta->GetMaximum() && fhGenGamEta->GetMaximum() >= fhGenEleEta->GetMaximum())
1736                 haxiseta->SetMaximum(fhGenGamEta->GetMaximum());
1737         else if(fhGenEtaEta->GetMaximum() >= fhGenPi0Eta->GetMaximum() && fhGenEtaEta->GetMaximum() >= fhGenGamEta->GetMaximum() &&
1738                         fhGenEtaEta->GetMaximum() >= fhGenOmegaEta->GetMaximum() && fhGenEtaEta->GetMaximum() >= fhGenEleEta->GetMaximum())
1739                 haxiseta->SetMaximum(fhGenEtaEta->GetMaximum());
1740         else if(fhGenOmegaEta->GetMaximum() >= fhGenPi0Eta->GetMaximum() && fhGenOmegaEta->GetMaximum() >= fhGenEtaEta->GetMaximum() &&
1741                         fhGenOmegaEta->GetMaximum() >= fhGenGamEta->GetMaximum() && fhGenOmegaEta->GetMaximum() >= fhGenEleEta->GetMaximum())
1742                 haxiseta->SetMaximum(fhGenOmegaEta->GetMaximum());
1743         else if(fhGenEleEta->GetMaximum() >= fhGenPi0Eta->GetMaximum() && fhGenEleEta->GetMaximum() >= fhGenEtaEta->GetMaximum() &&
1744                         fhGenEleEta->GetMaximum() >= fhGenOmegaEta->GetMaximum() && fhGenEleEta->GetMaximum() >= fhGenGamEta->GetMaximum())
1745                 haxiseta->SetMaximum(fhGenEleEta->GetMaximum());
1746
1747         haxiseta->Draw("axis");
1748         fhGenPi0Eta->Draw("same");
1749         fhGenGamEta->Draw("same");
1750         fhGenEtaEta->Draw("same");
1751         fhGenOmegaEta->Draw("same");
1752         fhGenEleEta->Draw("same");
1753
1754
1755         c10->cd(3) ;
1756         gPad->SetLogy();
1757         TH1F * haxisphi  = (TH1F*) fhGenPi0Phi->Clone("axisphi");
1758         haxisphi->SetTitle("Generated Particles #phi, |#eta| < 1");
1759         fhGenPi0Phi->SetLineColor(1);
1760         fhGenGamPhi->SetLineColor(4);
1761         fhGenEtaPhi->SetLineColor(2);
1762         fhGenOmegaPhi->SetLineColor(7);
1763         fhGenElePhi->SetLineColor(6);
1764         //Select the maximum of the histogram to show all lines.
1765         if(fhGenPi0Phi->GetMaximum() >= fhGenGamPhi->GetMaximum() && fhGenPi0Phi->GetMaximum() >= fhGenEtaPhi->GetMaximum() &&
1766            fhGenPi0Phi->GetMaximum() >= fhGenOmegaPhi->GetMaximum() && fhGenPi0Phi->GetMaximum() >= fhGenElePhi->GetMaximum())
1767                 haxisphi->SetMaximum(fhGenPi0Phi->GetMaximum());
1768         else if(fhGenGamPhi->GetMaximum() >= fhGenPi0Phi->GetMaximum() && fhGenGamPhi->GetMaximum() >= fhGenEtaPhi->GetMaximum() &&
1769                         fhGenGamPhi->GetMaximum() >= fhGenOmegaPhi->GetMaximum() && fhGenGamPhi->GetMaximum() >= fhGenElePhi->GetMaximum())
1770                 haxisphi->SetMaximum(fhGenGamPhi->GetMaximum());
1771         else if(fhGenEtaPhi->GetMaximum() >= fhGenPi0Phi->GetMaximum() && fhGenEtaPhi->GetMaximum() >= fhGenGamPhi->GetMaximum() &&
1772                         fhGenEtaPhi->GetMaximum() >= fhGenOmegaPhi->GetMaximum() && fhGenEtaPhi->GetMaximum() >= fhGenElePhi->GetMaximum())
1773                 haxisphi->SetMaximum(fhGenEtaPhi->GetMaximum());
1774         else if(fhGenOmegaPhi->GetMaximum() >= fhGenPi0Phi->GetMaximum() && fhGenOmegaPhi->GetMaximum() >= fhGenEtaPhi->GetMaximum() &&
1775                         fhGenOmegaPhi->GetMaximum() >= fhGenGamPhi->GetMaximum() && fhGenOmegaPhi->GetMaximum() >= fhGenElePhi->GetMaximum())
1776                 haxisphi->SetMaximum(fhGenOmegaPhi->GetMaximum());
1777         else if(fhGenElePhi->GetMaximum() >= fhGenPi0Phi->GetMaximum() && fhGenElePhi->GetMaximum() >= fhGenEtaPhi->GetMaximum() &&
1778                         fhGenElePhi->GetMaximum() >= fhGenOmegaPhi->GetMaximum() && fhGenElePhi->GetMaximum() >= fhGenGamPhi->GetMaximum())
1779                 haxisphi->SetMaximum(fhGenElePhi->GetMaximum());
1780
1781         haxisphi->Draw("axis");
1782         fhGenPi0Phi->Draw("same");
1783         fhGenGamPhi->Draw("same");
1784         fhGenEtaPhi->Draw("same");
1785         fhGenOmegaPhi->Draw("same");
1786         fhGenElePhi->Draw("same");
1787
1788         sprintf(name,"QA_%s_GeneratedDistributions.eps",fCalorimeter.Data());
1789         c10->Print(name);
1790
1791
1792         //Reconstructed clusters depending on its original particle.
1793         //printf("c1\n");
1794         sprintf(cname,"QA_%s_recgenid",fCalorimeter.Data());
1795         TCanvas  * c11 = new TCanvas(cname, "Reconstructed particles, function of their original particle ID", 400, 400) ;
1796         c11->Divide(2, 2);
1797
1798
1799         c11->cd(1) ;
1800         gPad->SetLogy();
1801         TH1F * hGamE   = (TH1F*) fhGamE->ProjectionY("hGamE",-1,-1);
1802         TH1F * hPi0E   = (TH1F*) fhPi0E->ProjectionY("hPi0E",-1,-1);
1803         TH1F * hEleE   = (TH1F*) fhEleE->ProjectionY("hEleE",-1,-1);
1804         TH1F * hNeHadE = (TH1F*) fhNeHadE->ProjectionY("hNeHadE",-1,-1);
1805         TH1F * hChHadE = (TH1F*) fhNeHadE->ProjectionY("hChHadE",-1,-1);
1806         TH1F * haxisE  = (TH1F*) hPi0E->Clone("axisE");
1807         haxisE->SetTitle("Reconstructed particles E, function of their original particle ID");
1808         hPi0E->SetLineColor(1);
1809         hGamE->SetLineColor(4);
1810         hNeHadE->SetLineColor(2);
1811         hChHadE->SetLineColor(7);
1812         hEleE->SetLineColor(6);
1813
1814         //Select the maximum of the histogram to show all lines.
1815         if(hPi0E->GetMaximum() >= hGamE->GetMaximum() && hPi0E->GetMaximum() >= hNeHadE->GetMaximum() &&
1816            hPi0E->GetMaximum() >= hChHadE->GetMaximum() && hPi0E->GetMaximum() >= hEleE->GetMaximum())
1817                 haxisE->SetMaximum(hPi0E->GetMaximum());
1818         else if(hGamE->GetMaximum() >= hPi0E->GetMaximum() && hGamE->GetMaximum() >= hNeHadE->GetMaximum() &&
1819                         hGamE->GetMaximum() >= hChHadE->GetMaximum() && hGamE->GetMaximum() >= hEleE->GetMaximum())
1820                 haxisE->SetMaximum(hGamE->GetMaximum());
1821         else if(hNeHadE->GetMaximum() >= hPi0E->GetMaximum() && hNeHadE->GetMaximum() >= hGamE->GetMaximum() &&
1822                         hNeHadE->GetMaximum() >= hChHadE->GetMaximum() && hNeHadE->GetMaximum() >= hEleE->GetMaximum())
1823                 haxisE->SetMaximum(hNeHadE->GetMaximum());
1824         else if(hChHadE->GetMaximum() >= hPi0E->GetMaximum() && hChHadE->GetMaximum() >= hNeHadE->GetMaximum() &&
1825                         hChHadE->GetMaximum() >= hGamE->GetMaximum() && hChHadE->GetMaximum() >= hEleE->GetMaximum())
1826                 haxisE->SetMaximum(hChHadE->GetMaximum());
1827         else if(hEleE->GetMaximum() >= hPi0E->GetMaximum() && hEleE->GetMaximum() >= hNeHadE->GetMaximum() &&
1828                         hEleE->GetMaximum() >= hChHadE->GetMaximum() && hEleE->GetMaximum() >= hGamE->GetMaximum())
1829                 haxisE->SetMaximum(hEleE->GetMaximum());
1830
1831         haxisE->Draw("axis");
1832         hPi0E->Draw("same");
1833         hGamE->Draw("same");
1834         hNeHadE->Draw("same");
1835         hChHadE->Draw("same");
1836         hEleE->Draw("same");
1837
1838         TLegend pLegend2(0.75,0.45,0.9,0.8);
1839         pLegend2.SetTextSize(0.06);
1840         pLegend2.AddEntry(hPi0E,"#pi^{0}","L");
1841         pLegend2.AddEntry(hGamE,"#gamma","L");
1842         pLegend2.AddEntry(hEleE,"e^{#pm}","L");
1843         pLegend2.AddEntry(hChHadE,"h^{#pm}","L");
1844         pLegend2.AddEntry(hNeHadE,"h^{0}","L");
1845         pLegend2.SetFillColor(10);
1846         pLegend2.SetBorderSize(1);
1847         pLegend2.Draw();
1848
1849
1850         c11->cd(2) ;
1851         gPad->SetLogy();
1852         //printf("%s, %s, %s, %s, %s\n",fhGamPt->GetName(),fhPi0Pt->GetName(),fhElePt->GetName(),fhNeHadPt->GetName(), fhChHadPt->GetName());
1853         TH1F * hGamPt   = (TH1F*) fhGamPt->ProjectionY("hGamPt",-1,-1);
1854         TH1F * hPi0Pt   = (TH1F*) fhPi0Pt->ProjectionY("hPi0Pt",-1,-1);
1855         TH1F * hElePt   = (TH1F*) fhElePt->ProjectionY("hElePt",-1,-1);
1856         TH1F * hNeHadPt = (TH1F*) fhNeHadPt->ProjectionY("hNeHadPt",-1,-1);
1857         TH1F * hChHadPt = (TH1F*) fhNeHadPt->ProjectionY("hChHadPt",-1,-1);
1858         haxispt  = (TH1F*) hPi0Pt->Clone("axispt");
1859         haxispt->SetTitle("Reconstructed particles p_{T}, function of their original particle ID");
1860         hPi0Pt->SetLineColor(1);
1861         hGamPt->SetLineColor(4);
1862         hNeHadPt->SetLineColor(2);
1863         hChHadPt->SetLineColor(7);
1864         hElePt->SetLineColor(6);
1865
1866         //Select the maximum of the histogram to show all lines.
1867         if(hPi0Pt->GetMaximum() >= hGamPt->GetMaximum() && hPi0Pt->GetMaximum() >= hNeHadPt->GetMaximum() &&
1868            hPi0Pt->GetMaximum() >= hChHadPt->GetMaximum() && hPi0Pt->GetMaximum() >= hElePt->GetMaximum())
1869                 haxispt->SetMaximum(hPi0Pt->GetMaximum());
1870         else if(hGamPt->GetMaximum() >= hPi0Pt->GetMaximum() && hGamPt->GetMaximum() >= hNeHadPt->GetMaximum() &&
1871                         hGamPt->GetMaximum() >= hChHadPt->GetMaximum() && hGamPt->GetMaximum() >= hElePt->GetMaximum())
1872                 haxispt->SetMaximum(hGamPt->GetMaximum());
1873         else if(hNeHadPt->GetMaximum() >= hPi0Pt->GetMaximum() && hNeHadPt->GetMaximum() >= hGamPt->GetMaximum() &&
1874                         hNeHadPt->GetMaximum() >= hChHadPt->GetMaximum() && hNeHadPt->GetMaximum() >= hElePt->GetMaximum())
1875                 haxispt->SetMaximum(hNeHadPt->GetMaximum());
1876         else if(hChHadPt->GetMaximum() >= hPi0Pt->GetMaximum() && hChHadPt->GetMaximum() >= hNeHadPt->GetMaximum() &&
1877                         hChHadPt->GetMaximum() >= hGamPt->GetMaximum() && hChHadPt->GetMaximum() >= hElePt->GetMaximum())
1878                 haxispt->SetMaximum(hChHadPt->GetMaximum());
1879         else if(hElePt->GetMaximum() >= hPi0Pt->GetMaximum() && hElePt->GetMaximum() >= hNeHadPt->GetMaximum() &&
1880                         hElePt->GetMaximum() >= hChHadPt->GetMaximum() && hElePt->GetMaximum() >= hGamPt->GetMaximum())
1881                 haxispt->SetMaximum(hElePt->GetMaximum());
1882
1883         haxispt->Draw("axis");
1884         hPi0Pt->Draw("same");
1885         hGamPt->Draw("same");
1886         hNeHadPt->Draw("same");
1887         hChHadPt->Draw("same");
1888         hElePt->Draw("same");
1889
1890
1891         c11->cd(3) ;
1892         gPad->SetLogy();
1893
1894         TH1F * hGamEta   = (TH1F*) fhGamEta->ProjectionY("hGamEta",-1,-1);
1895         TH1F * hPi0Eta   = (TH1F*) fhPi0Eta->ProjectionY("hPi0Eta",-1,-1);
1896         TH1F * hEleEta   = (TH1F*) fhEleEta->ProjectionY("hEleEta",-1,-1);
1897         TH1F * hNeHadEta = (TH1F*) fhNeHadEta->ProjectionY("hNeHadEta",-1,-1);
1898         TH1F * hChHadEta = (TH1F*) fhNeHadEta->ProjectionY("hChHadEta",-1,-1);
1899         haxiseta  = (TH1F*) hPi0Eta->Clone("axiseta");
1900         haxiseta->SetTitle("Reconstructed particles #eta, function of their original particle ID");
1901         hPi0Eta->SetLineColor(1);
1902         hGamEta->SetLineColor(4);
1903         hNeHadEta->SetLineColor(2);
1904         hChHadEta->SetLineColor(7);
1905         hEleEta->SetLineColor(6);
1906         //Select the maximum of the histogram to show all lines.
1907         if(hPi0Eta->GetMaximum() >= hGamEta->GetMaximum() && hPi0Eta->GetMaximum() >= hNeHadEta->GetMaximum() &&
1908            hPi0Eta->GetMaximum() >= hChHadEta->GetMaximum() && hPi0Eta->GetMaximum() >= hEleEta->GetMaximum())
1909                 haxiseta->SetMaximum(hPi0Eta->GetMaximum());
1910         else if(hGamEta->GetMaximum() >= hPi0Eta->GetMaximum() && hGamEta->GetMaximum() >= hNeHadEta->GetMaximum() &&
1911                         hGamEta->GetMaximum() >= hChHadEta->GetMaximum() && hGamEta->GetMaximum() >= hEleEta->GetMaximum())
1912                 haxiseta->SetMaximum(hGamEta->GetMaximum());
1913         else if(hNeHadEta->GetMaximum() >= hPi0Eta->GetMaximum() && hNeHadEta->GetMaximum() >= hGamEta->GetMaximum() &&
1914                         hNeHadEta->GetMaximum() >= hChHadEta->GetMaximum() && hNeHadEta->GetMaximum() >= hEleEta->GetMaximum())
1915                 haxiseta->SetMaximum(hNeHadEta->GetMaximum());
1916         else if(hChHadEta->GetMaximum() >= hPi0Eta->GetMaximum() && hChHadEta->GetMaximum() >= hNeHadEta->GetMaximum() &&
1917                         hChHadEta->GetMaximum() >= hGamEta->GetMaximum() && hChHadEta->GetMaximum() >= hEleEta->GetMaximum())
1918                 haxiseta->SetMaximum(hChHadEta->GetMaximum());
1919         else if(hEleEta->GetMaximum() >= hPi0Eta->GetMaximum() && hEleEta->GetMaximum() >= hNeHadEta->GetMaximum() &&
1920                         hEleEta->GetMaximum() >= hChHadEta->GetMaximum() && hEleEta->GetMaximum() >= hGamEta->GetMaximum())
1921                 haxiseta->SetMaximum(hEleEta->GetMaximum());
1922
1923         haxiseta->Draw("axis");
1924         hPi0Eta->Draw("same");
1925         hGamEta->Draw("same");
1926         hNeHadEta->Draw("same");
1927         hChHadEta->Draw("same");
1928         hEleEta->Draw("same");
1929
1930
1931         c11->cd(4) ;
1932         gPad->SetLogy();
1933         TH1F * hGamPhi   = (TH1F*) fhGamPhi->ProjectionY("hGamPhi",-1,-1);
1934         TH1F * hPi0Phi   = (TH1F*) fhPi0Phi->ProjectionY("hPi0Phi",-1,-1);
1935         TH1F * hElePhi   = (TH1F*) fhElePhi->ProjectionY("hElePhi",-1,-1);
1936         TH1F * hNeHadPhi = (TH1F*) fhNeHadPhi->ProjectionY("hNeHadPhi",-1,-1);
1937         TH1F * hChHadPhi = (TH1F*) fhNeHadPhi->ProjectionY("hChHadPhi",-1,-1);
1938         haxisphi  = (TH1F*) hPi0Phi->Clone("axisphi");
1939         haxisphi->SetTitle("Reconstructed particles #phi, function of their original particle ID");
1940
1941         hPi0Phi->SetLineColor(1);
1942         hGamPhi->SetLineColor(4);
1943         hNeHadPhi->SetLineColor(2);
1944         hChHadPhi->SetLineColor(7);
1945         hElePhi->SetLineColor(6);
1946         //Select the maximum of the histogram to show all lines.
1947         if(hPi0Phi->GetMaximum() >= hGamPhi->GetMaximum() && hPi0Phi->GetMaximum() >= hNeHadPhi->GetMaximum() &&
1948            hPi0Phi->GetMaximum() >= hChHadPhi->GetMaximum() && hPi0Phi->GetMaximum() >= hElePhi->GetMaximum())
1949                 haxisphi->SetMaximum(hPi0Phi->GetMaximum());
1950         else if(hGamPhi->GetMaximum() >= hPi0Phi->GetMaximum() && hGamPhi->GetMaximum() >= hNeHadPhi->GetMaximum() &&
1951                         hGamPhi->GetMaximum() >= hChHadPhi->GetMaximum() && hGamPhi->GetMaximum() >= hElePhi->GetMaximum())
1952                 haxisphi->SetMaximum(hGamPhi->GetMaximum());
1953         else if(hNeHadPhi->GetMaximum() >= hPi0Phi->GetMaximum() && hNeHadPhi->GetMaximum() >= hGamPhi->GetMaximum() &&
1954                         hNeHadPhi->GetMaximum() >= hChHadPhi->GetMaximum() && hNeHadPhi->GetMaximum() >= hElePhi->GetMaximum())
1955                 haxisphi->SetMaximum(hNeHadPhi->GetMaximum());
1956         else if(hChHadPhi->GetMaximum() >= hPi0Phi->GetMaximum() && hChHadPhi->GetMaximum() >= hNeHadPhi->GetMaximum() &&
1957                         hChHadPhi->GetMaximum() >= hGamPhi->GetMaximum() && hChHadPhi->GetMaximum() >= hElePhi->GetMaximum())
1958                 haxisphi->SetMaximum(hChHadPhi->GetMaximum());
1959         else if(hElePhi->GetMaximum() >= hPi0Phi->GetMaximum() && hElePhi->GetMaximum() >= hNeHadPhi->GetMaximum() &&
1960                         hElePhi->GetMaximum() >= hChHadPhi->GetMaximum() && hElePhi->GetMaximum() >= hGamPhi->GetMaximum())
1961                 haxisphi->SetMaximum(hElePhi->GetMaximum());
1962
1963         haxisphi->Draw("axis");
1964         hPi0Phi->Draw("same");
1965         hGamPhi->Draw("same");
1966         hNeHadPhi->Draw("same");
1967         hChHadPhi->Draw("same");
1968         hElePhi->Draw("same");
1969
1970         sprintf(name,"QA_%s_RecDistributionsGenID.eps",fCalorimeter.Data());
1971         c11->Print(name);
1972
1973
1974         //Ratio reconstructed clusters / generated particles in acceptance, for different particle ID
1975         //printf("c12\n");
1976
1977         TH1F *  hPi0EClone   = (TH1F*)   hPi0E  ->Clone("hPi0EClone");
1978         TH1F *  hGamEClone   = (TH1F*)   hGamE  ->Clone("hGamEClone");
1979         TH1F *  hPi0PtClone  = (TH1F*)   hPi0Pt ->Clone("hPi0PtClone");
1980         TH1F *  hGamPtClone  = (TH1F*)   hGamPt ->Clone("hGamPtClone");
1981         TH1F *  hPi0EtaClone = (TH1F*)   hPi0Eta->Clone("hPi0EtaClone");
1982         TH1F *  hGamEtaClone = (TH1F*)   hGamEta->Clone("hGamEtaClone");
1983         TH1F *  hPi0PhiClone = (TH1F*)   hPi0Phi->Clone("hPi0PhiClone");
1984         TH1F *  hGamPhiClone = (TH1F*)   hGamPhi->Clone("hGamPhiClone");
1985         sprintf(cname,"QA_%s_recgenidratio",fCalorimeter.Data());
1986         TCanvas  * c12 = new TCanvas(cname, "Ratio reconstructed clusters / generated particles in acceptance, for different particle ID", 400, 400) ;
1987         c12->Divide(2, 2);
1988
1989         c12->cd(1) ;
1990         gPad->SetLogy();
1991         haxisE->SetTitle("Ratio reconstructed clusters / generated particles in acceptance, for different particle ID");
1992         hPi0EClone->Divide(fhGenPi0AccE);
1993         hGamEClone->Divide(fhGenGamAccE);
1994         haxisE->SetMaximum(1.2);
1995         haxisE->SetYTitle("ratio = rec/gen");
1996         haxisE->Draw("axis");
1997         hPi0EClone->SetLineColor(1);
1998         hGamEClone->SetLineColor(4);
1999         hPi0EClone->Draw("same");
2000         hGamEClone->Draw("same");
2001
2002         TLegend pLegend3(0.75,0.45,0.9,0.8);
2003         pLegend3.SetTextSize(0.06);
2004         pLegend3.AddEntry(hPi0EClone,"#pi^{0}","L");
2005         pLegend3.AddEntry(hGamEClone,"#gamma","L");
2006         pLegend3.SetFillColor(10);
2007         pLegend3.SetBorderSize(1);
2008         pLegend3.Draw();
2009
2010         c12->cd(2) ;
2011         gPad->SetLogy();
2012         haxispt->SetTitle("Ratio reconstructed clusters / generated particles in acceptance, for different particle ID");
2013         hPi0PtClone->Divide(fhGenPi0AccPt);
2014         hGamPtClone->Divide(fhGenGamAccPt);
2015         haxispt->SetMaximum(1.2);
2016         haxispt->SetYTitle("ratio = rec/gen");
2017         haxispt->Draw("axis");
2018         hPi0PtClone->SetLineColor(1);
2019         hGamPtClone->SetLineColor(4);
2020         hPi0PtClone->Draw("same");
2021         hGamPtClone->Draw("same");
2022
2023         c12->cd(3) ;
2024         gPad->SetLogy();
2025
2026         haxiseta->SetTitle("Ratio reconstructed clusters / generated particles in acceptance, for different particle ID");
2027         hPi0EtaClone->Divide(fhGenPi0AccEta);
2028         hGamEtaClone->Divide(fhGenGamAccEta);
2029         haxiseta->SetMaximum(1.2);
2030         haxiseta->SetYTitle("ratio = rec/gen");
2031         haxiseta->Draw("axis");
2032         hPi0EtaClone->SetLineColor(1);
2033         hGamEtaClone->SetLineColor(4);
2034         hPi0EtaClone->Draw("same");
2035         hGamEtaClone->Draw("same");
2036
2037
2038         c12->cd(4) ;
2039         gPad->SetLogy();
2040         haxisphi->SetTitle("Ratio reconstructed clusters / generated particles in acceptance, for different particle ID");
2041         hPi0PhiClone->Divide(fhGenPi0AccPhi);
2042         hGamPhiClone->Divide(fhGenGamAccPhi);
2043         haxisphi->SetYTitle("ratio = rec/gen");
2044         haxisphi->SetMaximum(1.2);
2045         haxisphi->Draw("axis");
2046         hPi0PhiClone->SetLineColor(1);
2047         hGamPhiClone->SetLineColor(4);
2048         hPi0PhiClone->Draw("same");
2049         hGamPhiClone->Draw("same");
2050
2051
2052         sprintf(name,"QA_%s_EfficiencyGenID.eps",fCalorimeter.Data());
2053         c12->Print(name);
2054
2055
2056         //Reconstructed distributions
2057         //printf("c1\n");
2058         sprintf(cname,"QA_%s_vertex",fCalorimeter.Data());
2059         TCanvas  * c13 = new TCanvas(cname, "Particle vertex", 400, 400) ;
2060         c13->Divide(2, 2);
2061
2062         c13->cd(1) ;
2063         //gPad->SetLogy();
2064         fhEMVxyz->Draw();
2065
2066         c13->cd(2) ;
2067         //gPad->SetLogy();
2068         fhHaVxyz->Draw();
2069
2070         c13->cd(3) ;
2071         gPad->SetLogy();
2072         TH1F * hEMR = (TH1F*) fhEMR->ProjectionY("hEM",-1,-1);
2073         hEMR->SetLineColor(4);
2074         hEMR->Draw();
2075
2076         c13->cd(4) ;
2077         gPad->SetLogy();
2078         TH1F * hHaR = (TH1F*) fhHaR->ProjectionY("hHa",-1,-1);
2079         hHaR->SetLineColor(4);
2080         hHaR->Draw();
2081
2082
2083         sprintf(name,"QA_%s_ParticleVertex.eps",fCalorimeter.Data());
2084         c13->Print(name);
2085
2086
2087         //Track-matching distributions
2088         if(!strcmp(GetReader()->GetInputEvent()->GetName(),"AliESDEvent")){
2089                 sprintf(cname,"QA_%s_trkmatch",fCalorimeter.Data());
2090                 TCanvas *cme = new TCanvas(cname,"Track-matching distributions", 400, 400);
2091                 cme->Divide(2,2);
2092
2093                 cme->cd(1);
2094                 fh1pOverE->Draw();
2095
2096                 cme->cd(2);
2097                 fh1dR->Draw();
2098
2099                 cme->cd(3);
2100                 fh2MatchdEdx->Draw();
2101
2102                 cme->cd(4);
2103                 fh2EledEdx->Draw();
2104
2105                 sprintf(name,"QA_%s_TrackMatchingEleDist.eps",fCalorimeter.Data());
2106                 cme->Print(name);
2107         }
2108         char line[1024] ;
2109         sprintf(line, ".!tar -zcf QA_%s_%s.tar.gz *%s*.eps", fCalorimeter.Data(), GetName(),fCalorimeter.Data()) ;
2110         gROOT->ProcessLine(line);
2111         sprintf(line, ".!rm -fR *.eps");
2112         gROOT->ProcessLine(line);
2113
2114         printf("AliAnaCalorimeterQA::Terminate() - !! All the eps files are in QA_%s_%s.tar.gz !!!\n",  fCalorimeter.Data(), GetName());
2115
2116 }