PWGGA/CaloTrackCorrelations/AliAnaChargedParticles.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
  16 //_________________________________________________________________________
  17 //
  18 // Class for track selection and identification (not done now)
  19 // Tracks from the CTS are kept in the AOD.
  20 // Few histograms produced.
  21 //
  22 //-- Author: Gustavo Conesa (INFN-LNF)
  23 //_________________________________________________________________________
  24
  25
  26 // --- ROOT system ---
  27 #include "TParticle.h"
  28 #include "TH2F.h"
  29
  30 //---- AliRoot system ----
  31 #include "AliAnaChargedParticles.h"
  32 #include "AliCaloTrackReader.h"
  33 #include "AliAODPWG4Particle.h"
  34 #include "AliStack.h"
  35 #include "AliFiducialCut.h"
  36 #include "AliVTrack.h"
  37 #include "AliAODMCParticle.h"
  38 #include "AliAODTrack.h"
  39 #include "AliAODEvent.h"
  40 #include "AliESDEvent.h"
  41
  42 ClassImp(AliAnaChargedParticles)
  43
  44 //__________________________________________________
  45   AliAnaChargedParticles::AliAnaChargedParticles() :
  46     AliAnaCaloTrackCorrBaseClass(),
  47     fFillPileUpHistograms(0),   fFillTrackBCHistograms(0),
  48     fFillVertexBC0Histograms(0),
  49     //Histograms
  50     fhNtracks(0),      fhPt(0),            fhPtNoCut(0),
  51     fhPtCutDCA(0),     fhPtCutDCABCOK(0),
  52     fhPhiNeg(0),       fhEtaNeg(0),
  53     fhPhiPos(0),       fhEtaPos(0),
  54     fhEtaPhiPos(0),    fhEtaPhiNeg(0),
  55     fhPtVtxOutBC0(0),  fhEtaPhiVtxOutBC0(0),
  56     fhPtVtxInBC0(0),   fhEtaPhiVtxInBC0(0),
  57     fhPtSPDRefit(0),         fhPtNoSPDRefit(0),         fhPtNoSPDNoRefit(0),
  58     fhEtaPhiSPDRefitPt02(0), fhEtaPhiNoSPDRefitPt02(0), fhEtaPhiNoSPDNoRefitPt02(0),
  59     fhEtaPhiSPDRefitPt3(0),  fhEtaPhiNoSPDRefitPt3(0),  fhEtaPhiNoSPDNoRefitPt3(0),
  60     //MC
  61     fhPtPion(0),       fhPhiPion(0),         fhEtaPion(0),
  62     fhPtProton(0),     fhPhiProton(0),       fhEtaProton(0),
  63     fhPtElectron(0),   fhPhiElectron(0),     fhEtaElectron(0),
  64     fhPtKaon(0),       fhPhiKaon(0),         fhEtaKaon(0),
  65     fhPtUnknown(0),    fhPhiUnknown(0),      fhEtaUnknown(0),
  66     fhMCPt(0),         fhMCPhi(0),           fhMCEta(0),
  67     fhMCRecPt(0),
  68     //TOF
  69     fhTOFSignal(0),    fhTOFSignalPtCut(0),  fhTOFSignalBCOK(0),
  70     fhPtTOFSignal(0),  fhPtTOFSignalDCACut(0),
  71     fhPtTOFSignalVtxOutBC0(0), fhPtTOFSignalVtxInBC0(0),
  72     fhPtTOFStatus0(0), fhEtaPhiTOFStatus0(0),
  73     fhEtaPhiTOFBC0(0), fhEtaPhiTOFBCPlus(0), fhEtaPhiTOFBCMinus(0),
  74     fhEtaPhiTOFBC0PileUpSPD(0),
  75     fhEtaPhiTOFBCPlusPileUpSPD(0),
  76     fhEtaPhiTOFBCMinusPileUpSPD(0),
  77     fhProductionVertexBC(0)
  78 {
  79   //Default Ctor
  80
  81   for(Int_t i = 0; i < 7; i++)
  82   {
  83     fhPtPileUp         [i] = 0;
  84     fhPtTOFSignalPileUp[i] = 0;
  85     fhPtTOFSignalVtxOutBC0PileUp[i] = 0;
  86     fhPtTOFSignalVtxInBC0PileUp [i] = 0;
  87     fhProductionVertexBCPileUp  [i] = 0;
  88   }
  89
  90   for(Int_t i = 0; i < 3; i++)
  91   {
  92     fhPtDCA               [i] = 0 ;
  93
  94     fhPtDCASPDRefit       [i] = 0 ;
  95     fhPtDCANoSPDRefit     [i] = 0 ;
  96     fhPtDCANoSPDNoRefit   [i] = 0 ;
  97
  98     fhPtDCAPileUp         [i] = 0 ;
  99     fhPtDCATOFBC0         [i] = 0 ;
 100     fhPtDCATOFBCOut       [i] = 0 ;
 101     fhPtDCAPileUpTOFBC0   [i] = 0 ;
 102     fhPtDCANoTOFHit       [i] = 0 ;
 103     fhPtDCAPileUpNoTOFHit [i] = 0 ;
 104
 105     fhPtDCAVtxOutBC0      [i] = 0 ;
 106     fhPtDCAVtxInBC0       [i] = 0 ;
 107     fhPtDCAVtxOutBC0PileUp[i] = 0 ;
 108     fhPtDCAVtxInBC0PileUp [i] = 0 ;
 109     fhPtDCAVtxOutBC0NoTOFHit[i] = 0 ;
 110     fhPtDCAVtxInBC0NoTOFHit [i] = 0 ;
 111     fhPtDCAVtxOutBC0PileUpNoTOFHit[i] = 0 ;
 112     fhPtDCAVtxInBC0PileUpNoTOFHit [i] = 0 ;
 113   }
 114
 115   //Initialize parameters
 116   InitParameters();
 117
 118 }
 119
 120 //_______________________________________________________
 121 TList *  AliAnaChargedParticles::GetCreateOutputObjects()
 122 {
 123   // Create histograms to be saved in output file and
 124   // store them in fOutputContainer
 125
 126
 127   TList * outputContainer = new TList() ;
 128   outputContainer->SetName("ExampleHistos") ;
 129
 130   Int_t nptbins  = GetHistogramRanges()->GetHistoPtBins(); Int_t nphibins = GetHistogramRanges()->GetHistoPhiBins(); Int_t netabins = GetHistogramRanges()->GetHistoEtaBins();
 131   Float_t ptmax  = GetHistogramRanges()->GetHistoPtMax();  Float_t phimax = GetHistogramRanges()->GetHistoPhiMax();  Float_t etamax = GetHistogramRanges()->GetHistoEtaMax();
 132   Float_t ptmin  = GetHistogramRanges()->GetHistoPtMin();  Float_t phimin = GetHistogramRanges()->GetHistoPhiMin();  Float_t etamin = GetHistogramRanges()->GetHistoEtaMin();
 133
 134   fhNtracks  = new TH1F ("hNtracks","# of tracks", 1000,0,1000);
 135   fhNtracks->SetXTitle("# of tracks");
 136   outputContainer->Add(fhNtracks);
 137
 138   fhPt  = new TH1F ("hPt","p_T distribution", nptbins,ptmin,ptmax);
 139   fhPt->SetXTitle("p_{T} (GeV/c)");
 140   outputContainer->Add(fhPt);
 141
 142   fhPtNoCut  = new TH1F ("hPtNoCut","p_T distribution, raw tracks", nptbins,ptmin,ptmax);
 143   fhPtNoCut->SetXTitle("p_{T} (GeV/c)");
 144   outputContainer->Add(fhPtNoCut);
 145
 146   fhPtCutDCA  = new TH1F ("hPtCutDCA","p_T distribution, cut DCA", nptbins,ptmin,ptmax);
 147   fhPtCutDCA->SetXTitle("p_{T} (GeV/c)");
 148   outputContainer->Add(fhPtCutDCA);
 149
 150   if(fFillVertexBC0Histograms)
 151   {
 152     fhPtCutDCABCOK  = new TH1F ("hPtCutDCABCOK","p_T distribution, DCA cut, track BC=0 or -100", nptbins,ptmin,ptmax);
 153     fhPtCutDCABCOK->SetXTitle("p_{T} (GeV/c)");
 154     outputContainer->Add(fhPtCutDCABCOK);
 155   }
 156
 157   fhPhiNeg  = new TH2F ("hPhiNegative","#phi of negative charges distribution",
 158                         nptbins,ptmin,ptmax, nphibins,phimin,phimax);
 159   fhPhiNeg->SetYTitle("#phi (rad)");
 160   fhPhiNeg->SetXTitle("p_{T} (GeV/c)");
 161   outputContainer->Add(fhPhiNeg);
 162
 163   fhEtaNeg  = new TH2F ("hEtaNegative","#eta of negative charges distribution",
 164                         nptbins,ptmin,ptmax, netabins,etamin,etamax);
 165   fhEtaNeg->SetYTitle("#eta ");
 166   fhEtaNeg->SetXTitle("p_{T} (GeV/c)");
 167   outputContainer->Add(fhEtaNeg);
 168
 169   fhPhiPos  = new TH2F ("hPhiPositive","#phi of positive charges distribution",
 170                         nptbins,ptmin,ptmax, nphibins,phimin,phimax);
 171   fhPhiPos->SetYTitle("#phi (rad)");
 172   fhPhiPos->SetXTitle("p_{T} (GeV/c)");
 173   outputContainer->Add(fhPhiPos);
 174
 175   fhEtaPos  = new TH2F ("hEtaPositive","#eta of positive charges distribution",
 176                         nptbins,ptmin,ptmax, netabins,etamin,etamax);
 177   fhEtaPos->SetYTitle("#eta ");
 178   fhEtaPos->SetXTitle("p_{T} (GeV/c)");
 179   outputContainer->Add(fhEtaPos);
 180
 181   fhEtaPhiPos  = new TH2F ("hEtaPhiPositive","pt/eta/phi of positive charge",netabins,etamin,etamax, nphibins,phimin,phimax);
 182   fhEtaPhiPos->SetXTitle("#eta ");
 183   fhEtaPhiPos->SetYTitle("#phi (rad)");
 184   outputContainer->Add(fhEtaPhiPos);
 185
 186   fhEtaPhiNeg  = new TH2F ("hEtaPhiNegative","eta vs phi of negative charge",netabins,etamin,etamax, nphibins,phimin,phimax);
 187   fhEtaPhiNeg->SetXTitle("#eta ");
 188   fhEtaPhiNeg->SetYTitle("#phi (rad)");
 189   outputContainer->Add(fhEtaPhiNeg);
 190
 191   if(fFillVertexBC0Histograms)
 192   {
 193     fhPtVtxOutBC0  = new TH1F ("hPtVtxOutBC0","p_T distribution, vertex in BC=0", nptbins,ptmin,ptmax);
 194     fhPtVtxOutBC0->SetXTitle("p_{T} (GeV/c)");
 195     outputContainer->Add(fhPtVtxOutBC0);
 196
 197     fhEtaPhiVtxOutBC0  = new TH2F ("hEtaPhiVtxOutBC0","eta vs phi of all charges with vertex in BC=0",netabins,etamin,etamax, nphibins,phimin,phimax);
 198     fhEtaPhiVtxOutBC0->SetXTitle("#eta ");
 199     fhEtaPhiVtxOutBC0->SetYTitle("#phi (rad)");
 200     outputContainer->Add(fhEtaPhiVtxOutBC0);
 201
 202     fhPtVtxInBC0  = new TH1F ("hPtVtxInBC0","p_T distribution, vertex in BC=0", nptbins,ptmin,ptmax);
 203     fhPtVtxInBC0->SetXTitle("p_{T} (GeV/c)");
 204     outputContainer->Add(fhPtVtxInBC0);
 205
 206     fhEtaPhiVtxInBC0  = new TH2F ("hEtaPhiVtxInBC0","eta vs phi of all charges with vertex in BC=0",netabins,etamin,etamax, nphibins,phimin,phimax);
 207     fhEtaPhiVtxInBC0->SetXTitle("#eta ");
 208     fhEtaPhiVtxInBC0->SetYTitle("#phi (rad)");
 209     outputContainer->Add(fhEtaPhiVtxInBC0);
 210   }
 211
 212   fhPtSPDRefit  = new TH1F ("hPtSPDRefit","p_T distribution of tracks with SPD and ITS refit", nptbins,ptmin,ptmax);
 213   fhPtSPDRefit->SetXTitle("p_{T} (GeV/c)");
 214   outputContainer->Add(fhPtSPDRefit);
 215
 216   fhEtaPhiSPDRefitPt02  = new TH2F ("hEtaPhiSPDRefitPt02","eta vs phi of tracks with SPD and ITS refit, p_{T} < 2 GeV/c",netabins,etamin,etamax, nphibins,phimin,phimax);
 217   fhEtaPhiSPDRefitPt02->SetXTitle("#eta ");
 218   fhEtaPhiSPDRefitPt02->SetYTitle("#phi (rad)");
 219   outputContainer->Add(fhEtaPhiSPDRefitPt02);
 220
 221   fhEtaPhiSPDRefitPt3  = new TH2F ("hEtaPhiSPDRefitPt3","eta vs phi of tracks with SPD and ITS refit, p_{T} > 3 GeV/c",netabins,etamin,etamax, nphibins,phimin,phimax);
 222   fhEtaPhiSPDRefitPt3->SetXTitle("#eta ");
 223   fhEtaPhiSPDRefitPt3->SetYTitle("#phi (rad)");
 224   outputContainer->Add(fhEtaPhiSPDRefitPt3);
 225
 226   fhPtNoSPDRefit  = new TH1F ("hPtNoSPDRefit","p_T distribution of constrained tracks no SPD and with ITSRefit", nptbins,ptmin,ptmax);
 227   fhPtNoSPDRefit->SetXTitle("p_{T} (GeV/c)");
 228   outputContainer->Add(fhPtNoSPDRefit);
 229
 230   fhEtaPhiNoSPDRefitPt02  = new TH2F ("hEtaPhiNoSPDRefitPt02","eta vs phi of constrained tracks no SPD and with ITSRefit, p_{T} < 2 GeV/c",netabins,etamin,etamax, nphibins,phimin,phimax);
 231   fhEtaPhiNoSPDRefitPt02->SetXTitle("#eta ");
 232   fhEtaPhiNoSPDRefitPt02->SetYTitle("#phi (rad)");
 233   outputContainer->Add(fhEtaPhiNoSPDRefitPt02);
 234
 235   fhEtaPhiNoSPDRefitPt3  = new TH2F ("hEtaPhiNoSPDRefitPt3","eta vs phi of of constrained tracks no SPD and with ITSRefit, p_{T} > 3 GeV/c",netabins,etamin,etamax, nphibins,phimin,phimax);
 236   fhEtaPhiNoSPDRefitPt3->SetXTitle("#eta ");
 237   fhEtaPhiNoSPDRefitPt3->SetYTitle("#phi (rad)");
 238   outputContainer->Add(fhEtaPhiNoSPDRefitPt3);
 239
 240   fhPtNoSPDNoRefit  = new TH1F ("hPtNoSPDNoRefit","p_T distribution of constrained tracks with no SPD requierement and without ITSRefit", nptbins,ptmin,ptmax);
 241   fhPtNoSPDNoRefit->SetXTitle("p_{T} (GeV/c)");
 242   outputContainer->Add(fhPtNoSPDNoRefit);
 243
 244   fhEtaPhiNoSPDNoRefitPt02  = new TH2F ("hEtaPhiNoSPDNoRefitPt02","eta vs phi of constrained tracks with no SPD requierement and without ITSRefit, p_{T} < 2 GeV/c",netabins,etamin,etamax, nphibins,phimin,phimax);
 245   fhEtaPhiNoSPDNoRefitPt02->SetXTitle("#eta ");
 246   fhEtaPhiNoSPDNoRefitPt02->SetYTitle("#phi (rad)");
 247   outputContainer->Add(fhEtaPhiNoSPDNoRefitPt02);
 248
 249   fhEtaPhiNoSPDNoRefitPt3  = new TH2F ("hEtaPhiNoSPDNoRefitPt3","eta vs phi of constrained tracks with no SPD requierement and without ITSRefit, p_{T} > 3 GeV/c",netabins,etamin,etamax, nphibins,phimin,phimax);
 250   fhEtaPhiNoSPDNoRefitPt3->SetXTitle("#eta ");
 251   fhEtaPhiNoSPDNoRefitPt3->SetYTitle("#phi (rad)");
 252   outputContainer->Add(fhEtaPhiNoSPDNoRefitPt3);
 253
 254   if(fFillVertexBC0Histograms)
 255   {
 256     fhProductionVertexBC      = new TH1F("hProductionVertexBC", "tracks production vertex bunch crossing ", 41 , -20 , 20 ) ;
 257     fhProductionVertexBC->SetYTitle("# tracks");
 258     fhProductionVertexBC->SetXTitle("Bunch crossing");
 259     outputContainer->Add(fhProductionVertexBC);
 260   }
 261
 262   Int_t ntofbins = 1000;
 263   Int_t mintof = -500;
 264   Int_t maxtof =  500;
 265
 266   fhTOFSignal  = new TH1F ("hTOFSignal","TOF signal", ntofbins,mintof,maxtof);
 267   fhTOFSignal->SetXTitle("TOF signal (ns)");
 268   outputContainer->Add(fhTOFSignal);
 269
 270   if(fFillTrackBCHistograms)
 271   {
 272     fhTOFSignalBCOK  = new TH1F ("hTOFSignalBCOK","TOF signal", ntofbins,mintof,maxtof);
 273     fhTOFSignalBCOK->SetXTitle("TOF signal (ns)");
 274     outputContainer->Add(fhTOFSignalBCOK);
 275   }
 276
 277   fhTOFSignalPtCut  = new TH1F ("hTOFSignalPtCut","TOF signal", ntofbins,mintof,maxtof);
 278   fhTOFSignalPtCut->SetXTitle("TOF signal (ns)");
 279   outputContainer->Add(fhTOFSignalPtCut);
 280
 281   fhPtTOFSignal  = new TH2F ("hPtTOFSignal","TOF signal", nptbins,ptmin,ptmax,ntofbins,mintof,maxtof);
 282   fhPtTOFSignal->SetYTitle("TOF signal (ns)");
 283   fhPtTOFSignal->SetXTitle("p_{T} (GeV/c)");
 284   outputContainer->Add(fhPtTOFSignal);
 285
 286   fhPtTOFSignalDCACut  = new TH2F ("hPtTOFSignalDCACut","TOF signal after DCA cut", nptbins,ptmin,ptmax,ntofbins,mintof,maxtof);
 287   fhPtTOFSignalDCACut->SetYTitle("TOF signal (ns)");
 288   fhPtTOFSignalDCACut->SetXTitle("p_{T} (GeV/c)");
 289   outputContainer->Add(fhPtTOFSignalDCACut);
 290
 291   if(fFillVertexBC0Histograms)
 292   {
 293     fhPtTOFSignalVtxOutBC0  = new TH2F ("hPtTOFSignalVtxOutBC0","TOF signal, vtx BC!=0", nptbins,ptmin,ptmax,ntofbins,mintof,maxtof);
 294     fhPtTOFSignalVtxOutBC0->SetYTitle("TOF signal (ns)");
 295     fhPtTOFSignalVtxOutBC0->SetXTitle("p_{T} (GeV/c)");
 296     outputContainer->Add(fhPtTOFSignalVtxOutBC0);
 297
 298     fhPtTOFSignalVtxInBC0  = new TH2F ("hPtTOFSignalVtxInBC0","TOF signal, vtx BC=0", nptbins,ptmin,ptmax,ntofbins,mintof,maxtof);
 299     fhPtTOFSignalVtxInBC0->SetYTitle("TOF signal (ns)");
 300     fhPtTOFSignalVtxInBC0->SetXTitle("p_{T} (GeV/c)");
 301     outputContainer->Add(fhPtTOFSignalVtxInBC0);
 302   }
 303
 304   if(fFillPileUpHistograms)
 305   {
 306     TString pileUpName[] = {"SPD","EMCAL","SPDOrEMCAL","SPDAndEMCAL","SPDAndNotEMCAL","EMCALAndNotSPD","NotSPDAndNotEMCAL"} ;
 307
 308     for(Int_t i = 0 ; i < 7 ; i++)
 309     {
 310       fhPtPileUp[i]  = new TH1F(Form("hPtPileUp%s",pileUpName[i].Data()),
 311                                 Form("Track p_{T} distribution, %s Pile-Up event",pileUpName[i].Data()),
 312                                 nptbins,ptmin,ptmax);
 313       fhPtPileUp[i]->SetXTitle("p_{T} (GeV/c)");
 314       outputContainer->Add(fhPtPileUp[i]);
 315
 316       fhPtTOFSignalPileUp[i]  = new TH2F(Form("hPtTOFSignalPileUp%s",pileUpName[i].Data()),
 317                                          Form("Track TOF vs p_{T} distribution, %s Pile-Up event",pileUpName[i].Data()),
 318                                          nptbins,ptmin,ptmax,ntofbins,mintof,maxtof);
 319       fhPtTOFSignalPileUp[i]->SetXTitle("p_{T} (GeV/c)");
 320       fhPtTOFSignalPileUp[i]->SetXTitle("TOF signal (ns)");
 321       outputContainer->Add(fhPtTOFSignalPileUp[i]);
 322
 323       if(fFillVertexBC0Histograms)
 324       {
 325         fhPtTOFSignalVtxOutBC0PileUp[i]  = new TH2F(Form("hPtTOFSignalVtxOutBC0PileUp%s",pileUpName[i].Data()),
 326                                            Form("Track TOF vs p_{T} distribution, %s Pile-Up event, vtx BC!=0",pileUpName[i].Data()),
 327                                            nptbins,ptmin,ptmax,ntofbins,mintof,maxtof);
 328         fhPtTOFSignalVtxOutBC0PileUp[i]->SetXTitle("p_{T} (GeV/c)");
 329         fhPtTOFSignalVtxOutBC0PileUp[i]->SetXTitle("TOF signal (ns)");
 330         outputContainer->Add(fhPtTOFSignalVtxOutBC0PileUp[i]);
 331
 332         fhPtTOFSignalVtxInBC0PileUp[i]  = new TH2F(Form("hPtTOFSignalVtxInBC0PileUp%s",pileUpName[i].Data()),
 333                                                     Form("Track TOF vs p_{T} distribution, %s Pile-Up event, vtx BC=0",pileUpName[i].Data()),
 334                                                     nptbins,ptmin,ptmax,ntofbins,mintof,maxtof);
 335         fhPtTOFSignalVtxInBC0PileUp[i]->SetXTitle("p_{T} (GeV/c)");
 336         fhPtTOFSignalVtxInBC0PileUp[i]->SetXTitle("TOF signal (ns)");
 337         outputContainer->Add(fhPtTOFSignalVtxInBC0PileUp[i]);
 338       }
 339
 340       if(fFillVertexBC0Histograms)
 341       {
 342         fhProductionVertexBCPileUp[i]      = new TH1F(Form("hProductionVertexBCPileUp%s",pileUpName[i].Data()),
 343                                                 Form("tracks production vertex bunch crossing, %s Pile-Up event",pileUpName[i].Data()),
 344                                                 41 , -20 , 20 ) ;
 345         fhProductionVertexBCPileUp[i]->SetYTitle("# tracks");
 346         fhProductionVertexBCPileUp[i]->SetXTitle("Bunch crossing");
 347         outputContainer->Add(fhProductionVertexBCPileUp[i]);
 348       }
 349     }
 350
 351
 352     if(fFillTrackBCHistograms)
 353     {
 354       fhEtaPhiTOFBC0  = new TH2F ("hEtaPhiTOFBC0","eta-phi for tracks with hit on TOF, and tof corresponding to BC=0",netabins,etamin,etamax, nphibins,phimin,phimax);
 355       fhEtaPhiTOFBC0->SetXTitle("#eta ");
 356       fhEtaPhiTOFBC0->SetYTitle("#phi (rad)");
 357       outputContainer->Add(fhEtaPhiTOFBC0);
 358
 359       fhEtaPhiTOFBCPlus  = new TH2F ("hEtaPhiTOFBCPlus","eta-phi for tracks with hit on TOF, and tof corresponding to BC>0",netabins,etamin,etamax, nphibins,phimin,phimax);
 360       fhEtaPhiTOFBCPlus->SetXTitle("#eta ");
 361       fhEtaPhiTOFBCPlus->SetYTitle("#phi (rad)");
 362       outputContainer->Add(fhEtaPhiTOFBCPlus);
 363
 364       fhEtaPhiTOFBCMinus  = new TH2F ("hEtaPhiTOFBCMinus","eta-phi for tracks with hit on TOF, and tof corresponding to BC<0",netabins,etamin,etamax, nphibins,phimin,phimax);
 365       fhEtaPhiTOFBCMinus->SetXTitle("#eta ");
 366       fhEtaPhiTOFBCMinus->SetYTitle("#phi (rad)");
 367       outputContainer->Add(fhEtaPhiTOFBCMinus);
 368
 369       if(fFillPileUpHistograms)
 370       {
 371         fhEtaPhiTOFBC0PileUpSPD  = new TH2F ("hEtaPhiTOFBC0PileUpSPD","eta-phi for tracks with hit on TOF, and tof corresponding to BC=0, SPD pile-up",netabins,etamin,etamax, nphibins,phimin,phimax);
 372         fhEtaPhiTOFBC0PileUpSPD->SetXTitle("#eta ");
 373         fhEtaPhiTOFBC0PileUpSPD->SetYTitle("#phi (rad)");
 374         outputContainer->Add(fhEtaPhiTOFBC0PileUpSPD);
 375
 376         fhEtaPhiTOFBCPlusPileUpSPD  = new TH2F ("hEtaPhiTOFBCPlusPileUpSPD","eta-phi for tracks with hit on TOF, and tof corresponding to BC>0, SPD pile-up",netabins,etamin,etamax, nphibins,phimin,phimax);
 377         fhEtaPhiTOFBCPlusPileUpSPD->SetXTitle("#eta ");
 378         fhEtaPhiTOFBCPlusPileUpSPD->SetYTitle("#phi (rad)");
 379         outputContainer->Add(fhEtaPhiTOFBCPlusPileUpSPD);
 380
 381         fhEtaPhiTOFBCMinusPileUpSPD  = new TH2F ("hEtaPhiTOFBCMinusPileUpSPD","eta-phi for tracks with hit on TOF, and tof corresponding to BC<0, SPD pile-up",netabins,etamin,etamax, nphibins,phimin,phimax);
 382         fhEtaPhiTOFBCMinusPileUpSPD->SetXTitle("#eta ");
 383         fhEtaPhiTOFBCMinusPileUpSPD->SetYTitle("#phi (rad)");
 384         outputContainer->Add(fhEtaPhiTOFBCMinusPileUpSPD);
 385       }
 386     }
 387
 388   }
 389
 390   fhPtTOFStatus0  = new TH1F ("hPtTOFStatus0","p_T distribution of tracks not hitting TOF", nptbins,ptmin,ptmax);
 391   fhPtTOFStatus0->SetXTitle("p_{T} (GeV/c)");
 392   outputContainer->Add(fhPtTOFStatus0);
 393
 394
 395   fhEtaPhiTOFStatus0  = new TH2F ("hEtaPhiTOFStatus0","eta-phi for tracks without hit on TOF",netabins,etamin,etamax, nphibins,phimin,phimax);
 396   fhEtaPhiTOFStatus0->SetXTitle("#eta ");
 397   fhEtaPhiTOFStatus0->SetYTitle("#phi (rad)");
 398   outputContainer->Add(fhEtaPhiTOFStatus0);
 399
 400   TString dcaName[] = {"xy","z","Cons"} ;
 401   Int_t ndcabins = 800;
 402   Int_t mindca = -4;
 403   Int_t maxdca =  4;
 404
 405   for(Int_t i = 0 ; i < 3 ; i++)
 406   {
 407
 408     fhPtDCA[i]  = new TH2F(Form("hPtDCA%s",dcaName[i].Data()),
 409                            Form("Track DCA%s vs p_{T} distribution",dcaName[i].Data()),
 410                            nptbins,ptmin,ptmax,ndcabins,mindca,maxdca);
 411     fhPtDCA[i]->SetXTitle("p_{T} (GeV/c)");
 412     fhPtDCA[i]->SetYTitle(Form("DCA_{%s}",dcaName[i].Data()));
 413     outputContainer->Add(fhPtDCA[i]);
 414
 415     fhPtDCASPDRefit[i]  = new TH2F(Form("hPtDCA%sSPDRefit",dcaName[i].Data()),
 416                                         Form("Track DCA%s vs p_{T} distribution of tracks with SPD and ITS refit",dcaName[i].Data()),
 417                                         nptbins,ptmin,ptmax,ndcabins,mindca,maxdca);
 418     fhPtDCASPDRefit[i]->SetXTitle("p_{T} (GeV/c)");
 419     fhPtDCASPDRefit[i]->SetYTitle(Form("DCA_{%s}",dcaName[i].Data()));
 420     outputContainer->Add(fhPtDCASPDRefit[i]);
 421
 422     fhPtDCANoSPDRefit[i]  = new TH2F(Form("hPtDCA%sNoSPDRefit",dcaName[i].Data()),
 423                                  Form("Track DCA%s vs p_{T} distributionof constrained tracks no SPD and with ITSRefit",dcaName[i].Data()),
 424                                  nptbins,ptmin,ptmax,ndcabins,mindca,maxdca);
 425     fhPtDCANoSPDRefit[i]->SetXTitle("p_{T} (GeV/c)");
 426     fhPtDCANoSPDRefit[i]->SetYTitle(Form("DCA_{%s}",dcaName[i].Data()));
 427     outputContainer->Add(fhPtDCANoSPDRefit[i]);
 428
 429     fhPtDCANoSPDNoRefit[i]  = new TH2F(Form("hPtDCA%sNoSPDNoRefit",dcaName[i].Data()),
 430                            Form("Track DCA%s vs p_{T} distribution, constrained tracks with no SPD requierement and without ITSRefit",dcaName[i].Data()),
 431                            nptbins,ptmin,ptmax,ndcabins,mindca,maxdca);
 432     fhPtDCANoSPDNoRefit[i]->SetXTitle("p_{T} (GeV/c)");
 433     fhPtDCANoSPDNoRefit[i]->SetYTitle(Form("DCA_{%s}",dcaName[i].Data()));
 434     outputContainer->Add(fhPtDCANoSPDNoRefit[i]);
 435
 436     if(fFillTrackBCHistograms)
 437     {
 438       fhPtDCATOFBC0[i]  = new TH2F(Form("hPtDCA%sTOFBC0",dcaName[i].Data()),
 439                                    Form("Track DCA%s vs p_{T} distribution, BC=0",dcaName[i].Data()),
 440                                    nptbins,ptmin,ptmax,ndcabins,mindca,maxdca);
 441       fhPtDCATOFBC0[i]->SetXTitle("p_{T} (GeV/c)");
 442       fhPtDCATOFBC0[i]->SetYTitle(Form("DCA_{%s}",dcaName[i].Data()));
 443       outputContainer->Add(fhPtDCATOFBC0[i]);
 444
 445       fhPtDCATOFBCOut[i]  = new TH2F(Form("hPtDCA%sTOFBCOut",dcaName[i].Data()),
 446                                      Form("Track DCA%s vs p_{T} distribution, BC!=0",dcaName[i].Data()),
 447                                      nptbins,ptmin,ptmax,ndcabins,mindca,maxdca);
 448       fhPtDCATOFBCOut[i]->SetXTitle("p_{T} (GeV/c)");
 449       fhPtDCATOFBCOut[i]->SetYTitle(Form("DCA_{%s}",dcaName[i].Data()));
 450       outputContainer->Add(fhPtDCATOFBCOut[i]);
 451     }
 452
 453     fhPtDCANoTOFHit[i]  = new TH2F(Form("hPtDCA%sNoTOFHit",dcaName[i].Data()),
 454                            Form("Track (no TOF hit) DCA%s vs p_{T} distribution",dcaName[i].Data()),
 455                            nptbins,ptmin,ptmax,ndcabins,mindca,maxdca);
 456     fhPtDCANoTOFHit[i]->SetXTitle("p_{T} (GeV/c)");
 457     fhPtDCANoTOFHit[i]->SetYTitle(Form("DCA_{%s}",dcaName[i].Data()));
 458     outputContainer->Add(fhPtDCANoTOFHit[i]);
 459
 460     if(fFillVertexBC0Histograms)
 461     {
 462       fhPtDCAVtxOutBC0[i]  = new TH2F(Form("hPtDCA%sVtxOutBC0",dcaName[i].Data()),
 463                                       Form("Track DCA%s vs p_{T} distribution, vertex with BC!=0",dcaName[i].Data()),
 464                                       nptbins,ptmin,ptmax,ndcabins,mindca,maxdca);
 465       fhPtDCAVtxOutBC0[i]->SetXTitle("p_{T} (GeV/c)");
 466       fhPtDCAVtxOutBC0[i]->SetYTitle(Form("DCA_{%s}",dcaName[i].Data()));
 467       outputContainer->Add(fhPtDCAVtxOutBC0[i]);
 468
 469       fhPtDCAVtxOutBC0NoTOFHit[i]  = new TH2F(Form("hPtDCA%sVtxOutBC0NoTOFHit",dcaName[i].Data()),
 470                                               Form("Track (no TOF hit) DCA%s vs p_{T} distribution, vertex with BC!=0",dcaName[i].Data()),
 471                                               nptbins,ptmin,ptmax,ndcabins,mindca,maxdca);
 472       fhPtDCAVtxOutBC0NoTOFHit[i]->SetXTitle("p_{T} (GeV/c)");
 473       fhPtDCAVtxOutBC0NoTOFHit[i]->SetYTitle(Form("DCA_{%s}",dcaName[i].Data()));
 474       outputContainer->Add(fhPtDCAVtxOutBC0NoTOFHit[i]);
 475
 476       fhPtDCAVtxInBC0[i]  = new TH2F(Form("hPtDCA%sVtxInBC0",dcaName[i].Data()),
 477                                       Form("Track DCA%s vs p_{T} distribution, vertex with BC==0",dcaName[i].Data()),
 478                                       nptbins,ptmin,ptmax,ndcabins,mindca,maxdca);
 479       fhPtDCAVtxInBC0[i]->SetXTitle("p_{T} (GeV/c)");
 480       fhPtDCAVtxInBC0[i]->SetYTitle(Form("DCA_{%s}",dcaName[i].Data()));
 481       outputContainer->Add(fhPtDCAVtxInBC0[i]);
 482
 483       fhPtDCAVtxInBC0NoTOFHit[i]  = new TH2F(Form("hPtDCA%sVtxInBC0NoTOFHit",dcaName[i].Data()),
 484                                               Form("Track (no TOF hit) DCA%s vs p_{T} distribution, vertex with BC==0",dcaName[i].Data()),
 485                                               nptbins,ptmin,ptmax,ndcabins,mindca,maxdca);
 486       fhPtDCAVtxInBC0NoTOFHit[i]->SetXTitle("p_{T} (GeV/c)");
 487       fhPtDCAVtxInBC0NoTOFHit[i]->SetYTitle(Form("DCA_{%s}",dcaName[i].Data()));
 488       outputContainer->Add(fhPtDCAVtxInBC0NoTOFHit[i]);
 489     }
 490
 491     if(fFillPileUpHistograms)
 492     {
 493       fhPtDCAPileUp[i]  = new TH2F(Form("hPtDCA%sPileUp",dcaName[i].Data()),
 494                              Form("Track DCA%s vs p_{T} distribution, SPD Pile-Up",dcaName[i].Data()),
 495                              nptbins,ptmin,ptmax,ndcabins,mindca,maxdca);
 496       fhPtDCAPileUp[i]->SetXTitle("p_{T} (GeV/c)");
 497       fhPtDCAPileUp[i]->SetYTitle(Form("DCA_{%s}",dcaName[i].Data()));
 498       outputContainer->Add(fhPtDCAPileUp[i]);
 499
 500       if(fFillTrackBCHistograms)
 501       {
 502         fhPtDCAPileUpTOFBC0[i]  = new TH2F(Form("hPtDCA%sPileUpTOFBC0",dcaName[i].Data()),
 503                                            Form("Track DCA%s vs p_{T} distribution",dcaName[i].Data()),
 504                                            nptbins,ptmin,ptmax,ndcabins,mindca,maxdca);
 505         fhPtDCAPileUpTOFBC0[i]->SetXTitle("p_{T} (GeV/c)");
 506         fhPtDCAPileUpTOFBC0[i]->SetYTitle(Form("DCA_{%s}",dcaName[i].Data()));
 507         outputContainer->Add(fhPtDCAPileUpTOFBC0[i]);
 508       }
 509
 510       fhPtDCAPileUpNoTOFHit[i]  = new TH2F(Form("hPtDCA%sPileUpNoTOFHit",dcaName[i].Data()),
 511                                      Form("Track (no TOF hit) DCA%s vs p_{T} distribution, SPD Pile-Up, vertex with BC!=0",dcaName[i].Data()),
 512                                      nptbins,ptmin,ptmax,ndcabins,mindca,maxdca);
 513       fhPtDCAPileUpNoTOFHit[i]->SetXTitle("p_{T} (GeV/c)");
 514       fhPtDCAPileUpNoTOFHit[i]->SetYTitle(Form("DCA_{%s}",dcaName[i].Data()));
 515       outputContainer->Add(fhPtDCAPileUpNoTOFHit[i]);
 516
 517       if(fFillVertexBC0Histograms)
 518       {
 519         fhPtDCAVtxOutBC0PileUp[i]  = new TH2F(Form("hPtDCA%sPileUpVtxOutBC0",dcaName[i].Data()),
 520                                               Form("Track DCA%s vs p_{T} distribution, SPD Pile-Up, vertex with BC!=0",dcaName[i].Data()),
 521                                               nptbins,ptmin,ptmax,ndcabins,mindca,maxdca);
 522         fhPtDCAVtxOutBC0PileUp[i]->SetXTitle("p_{T} (GeV/c)");
 523         fhPtDCAVtxOutBC0PileUp[i]->SetYTitle(Form("DCA_{%s}",dcaName[i].Data()));
 524         outputContainer->Add(fhPtDCAVtxOutBC0PileUp[i]);
 525
 526         fhPtDCAVtxOutBC0PileUpNoTOFHit[i]  = new TH2F(Form("hPtDCA%sVtxOutBC0PileUpNoTOFHit",dcaName[i].Data()),
 527                                                       Form("Track (no TOF hit) DCA%s vs p_{T} distribution, SPD Pile-Up, vertex with BC!=0",dcaName[i].Data()),
 528                                                       nptbins,ptmin,ptmax,ndcabins,mindca,maxdca);
 529         fhPtDCAVtxOutBC0PileUpNoTOFHit[i]->SetXTitle("p_{T} (GeV/c)");
 530         fhPtDCAVtxOutBC0PileUpNoTOFHit[i]->SetYTitle(Form("DCA_{%s}",dcaName[i].Data()));
 531         outputContainer->Add(fhPtDCAVtxOutBC0PileUpNoTOFHit[i]);
 532
 533         fhPtDCAVtxInBC0PileUp[i]  = new TH2F(Form("hPtDCA%sPileUpVtxInBC0",dcaName[i].Data()),
 534                                               Form("Track DCA%s vs p_{T} distribution, SPD Pile-Up,vertex with BC==0",dcaName[i].Data()),
 535                                               nptbins,ptmin,ptmax,ndcabins,mindca,maxdca);
 536         fhPtDCAVtxInBC0PileUp[i]->SetXTitle("p_{T} (GeV/c)");
 537         fhPtDCAVtxInBC0PileUp[i]->SetYTitle(Form("DCA_{%s}",dcaName[i].Data()));
 538         outputContainer->Add(fhPtDCAVtxInBC0PileUp[i]);
 539
 540         fhPtDCAVtxInBC0PileUpNoTOFHit[i]  = new TH2F(Form("hPtDCA%sVtxInBC0PileUpNoTOFHit",dcaName[i].Data()),
 541                                                       Form("Track (no TOF hit) DCA%s vs p_{T} distribution, SPD Pile-Up, vertex with BC==0",dcaName[i].Data()),
 542                                                       nptbins,ptmin,ptmax,ndcabins,mindca,maxdca);
 543         fhPtDCAVtxInBC0PileUpNoTOFHit[i]->SetXTitle("p_{T} (GeV/c)");
 544         fhPtDCAVtxInBC0PileUpNoTOFHit[i]->SetYTitle(Form("DCA_{%s}",dcaName[i].Data()));
 545         outputContainer->Add(fhPtDCAVtxInBC0PileUpNoTOFHit[i]);
 546
 547       }
 548     }
 549   }
 550
 551   if(IsDataMC())
 552   {
 553     fhPtPion  = new TH1F ("hPtMCPion","p_T distribution from #pi", nptbins,ptmin,ptmax);
 554     fhPtPion->SetXTitle("p_{T} (GeV/c)");
 555     outputContainer->Add(fhPtPion);
 556
 557     fhPhiPion  = new TH2F ("hPhiMCPion","#phi distribution from #pi",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
 558     fhPhiPion->SetXTitle("#phi (rad)");
 559     outputContainer->Add(fhPhiPion);
 560
 561     fhEtaPion  = new TH2F ("hEtaMCPion","#eta distribution from #pi",nptbins,ptmin,ptmax, netabins,etamin,etamax);
 562     fhEtaPion->SetXTitle("#eta ");
 563     outputContainer->Add(fhEtaPion);
 564
 565     fhPtProton  = new TH1F ("hPtMCProton","p_T distribution from proton", nptbins,ptmin,ptmax);
 566     fhPtProton->SetXTitle("p_{T} (GeV/c)");
 567     outputContainer->Add(fhPtProton);
 568
 569     fhPhiProton  = new TH2F ("hPhiMCProton","#phi distribution from proton",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
 570     fhPhiProton->SetXTitle("#phi (rad)");
 571     outputContainer->Add(fhPhiProton);
 572
 573     fhEtaProton  = new TH2F ("hEtaMCProton","#eta distribution from proton",nptbins,ptmin,ptmax, netabins,etamin,etamax);
 574     fhEtaProton->SetXTitle("#eta ");
 575     outputContainer->Add(fhEtaProton);
 576
 577     fhPtKaon  = new TH1F ("hPtMCKaon","p_T distribution from kaon", nptbins,ptmin,ptmax);
 578     fhPtKaon->SetXTitle("p_{T} (GeV/c)");
 579     outputContainer->Add(fhPtKaon);
 580
 581     fhPhiKaon  = new TH2F ("hPhiMCKaon","#phi distribution from kaon",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
 582     fhPhiKaon->SetXTitle("#phi (rad)");
 583     outputContainer->Add(fhPhiKaon);
 584
 585     fhEtaKaon  = new TH2F ("hEtaMCKaon","#eta distribution from kaon",nptbins,ptmin,ptmax, netabins,etamin,etamax);
 586     fhEtaKaon->SetXTitle("#eta ");
 587     outputContainer->Add(fhEtaKaon);
 588
 589     fhPtElectron  = new TH1F ("hPtMCElectron","p_T distribution from electron", nptbins,ptmin,ptmax);
 590     fhPtElectron->SetXTitle("p_{T} (GeV/c)");
 591     outputContainer->Add(fhPtElectron);
 592
 593     fhPhiElectron  = new TH2F ("hPhiMCElectron","#phi distribution from electron",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
 594     fhPhiElectron->SetXTitle("#phi (rad)");
 595     outputContainer->Add(fhPhiElectron);
 596
 597     fhEtaElectron  = new TH2F ("hEtaMCElectron","#eta distribution from electron",nptbins,ptmin,ptmax, netabins,etamin,etamax);
 598     fhEtaElectron->SetXTitle("#eta ");
 599     outputContainer->Add(fhEtaElectron);
 600
 601     fhPtUnknown  = new TH1F ("hPtMCUnknown","p_T distribution from unknown", nptbins,ptmin,ptmax);
 602     fhPtUnknown->SetXTitle("p_{T} (GeV/c)");
 603     outputContainer->Add(fhPtUnknown);
 604
 605     fhPhiUnknown  = new TH2F ("hPhiMCUnknown","#phi distribution from unknown",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
 606     fhPhiUnknown->SetXTitle("#phi (rad)");
 607     outputContainer->Add(fhPhiUnknown);
 608
 609     fhEtaUnknown  = new TH2F ("hEtaMCUnknown","#eta distribution from unknown",nptbins,ptmin,ptmax, netabins,etamin,etamax);
 610     fhEtaUnknown->SetXTitle("#eta ");
 611     outputContainer->Add(fhEtaUnknown);
 612
 613     fhMCPt  = new TH1F ("hMCPt","p_T distribution from MC", nptbins,ptmin,ptmax);
 614     fhMCPt->SetXTitle("p_{T} (GeV/c)");
 615     outputContainer->Add(fhMCPt);
 616
 617     fhMCPhi  = new TH2F ("hMCPhi","#phi distribution from MC",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
 618     fhMCPhi->SetXTitle("#phi (rad)");
 619     outputContainer->Add(fhMCPhi);
 620
 621     fhMCEta  = new TH2F ("hMCEta","#eta distribution from MC",nptbins,ptmin,ptmax,netabins,etamin,etamax);
 622     fhMCEta->SetXTitle("#eta (rad)");
 623     outputContainer->Add(fhMCEta);
 624
 625     fhMCRecPt  = new TH1F ("hMCRecPt","p_T distribution from Rec MC", nptbins,ptmin,ptmax);
 626     fhMCRecPt->SetXTitle("p_{T} (GeV/c)");
 627     outputContainer->Add(fhMCRecPt);
 628   }
 629
 630   return outputContainer;
 631
 632 }
 633
 634
 635 //___________________________________________
 636 void AliAnaChargedParticles::InitParameters()
 637 {
 638   //Initialize the parameters of the analysis.
 639   SetOutputAODClassName("AliAODPWG4Particle");
 640   SetOutputAODName("PWG4Particle");
 641
 642   AddToHistogramsName("AnaCharged_");
 643
 644
 645 }
 646
 647 //____________________________________________________________
 648 void AliAnaChargedParticles::Print(const Option_t * opt) const
 649 {
 650   //Print some relevant parameters set for the analysis
 651   if(! opt)
 652     return;
 653
 654   printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
 655   AliAnaCaloTrackCorrBaseClass::Print(" ");
 656
 657   printf("Min Pt = %3.2f\n", GetMinPt());
 658   printf("Max Pt = %3.2f\n", GetMaxPt());
 659
 660 }
 661
 662 //_________________________________
 663 void AliAnaChargedParticles::Init()
 664 {
 665   //Init
 666   //Do some checks
 667   if(!GetReader()->IsCTSSwitchedOn()){
 668     printf("AliAnaChargedParticles::Init() - STOP!: You want to use CTS tracks in analysis but not read!! \n!!Check the configuration file!!\n");
 669     abort();
 670   }
 671
 672 }
 673
 674 //_________________________________________________
 675 void  AliAnaChargedParticles::MakeAnalysisFillAOD()
 676 {
 677   //Do analysis and fill aods
 678   if(!GetCTSTracks() || GetCTSTracks()->GetEntriesFast() == 0) return ;
 679
 680   Int_t ntracks = GetCTSTracks()->GetEntriesFast();
 681   Double_t vert[3] = {0,0,0}; //vertex ;
 682
 683   //Some prints
 684   if(GetDebug() > 0)
 685     printf("AliAnaChargedParticles::MakeAnalysisFillAOD() - In CTS aod entries %d\n", ntracks);
 686
 687   AliVEvent * event = GetReader()->GetInputEvent();
 688
 689   Int_t vtxBC = GetReader()->GetVertexBC();
 690   if(!GetReader()->IsDCACutOn()) vtxBC = GetReader()->GetVertexBC(event->GetPrimaryVertex());
 691
 692   if(fFillVertexBC0Histograms)
 693   {
 694     fhProductionVertexBC->Fill(vtxBC);
 695     if(fFillPileUpHistograms)
 696     {
 697       if(GetReader()->IsPileUpFromSPD())               fhProductionVertexBCPileUp[0]->Fill(vtxBC);
 698       if(GetReader()->IsPileUpFromEMCal())             fhProductionVertexBCPileUp[1]->Fill(vtxBC);
 699       if(GetReader()->IsPileUpFromSPDOrEMCal())        fhProductionVertexBCPileUp[2]->Fill(vtxBC);
 700       if(GetReader()->IsPileUpFromSPDAndEMCal())       fhProductionVertexBCPileUp[3]->Fill(vtxBC);
 701       if(GetReader()->IsPileUpFromSPDAndNotEMCal())    fhProductionVertexBCPileUp[4]->Fill(vtxBC);
 702       if(GetReader()->IsPileUpFromEMCalAndNotSPD())    fhProductionVertexBCPileUp[5]->Fill(vtxBC);
 703       if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) fhProductionVertexBCPileUp[6]->Fill(vtxBC);
 704     }
 705   }
 706
 707   //printf("AliAnaChargedParticles::MakeAnalysisFillAOD() - primary vertex BC %d\n",vtxBC);
 708
 709   Double_t bz = event->GetMagneticField();
 710
 711   //Fill AODParticle with CTS aods
 712   Float_t pt  = 0;
 713   Float_t phi = 0;
 714   Float_t eta = 0;
 715   Int_t evtIndex = 0;
 716   for(Int_t i = 0; i < ntracks; i++)
 717   {
 718     AliVTrack * track =  (AliVTrack*) (GetCTSTracks()->At(i));
 719
 720     pt  = track->Pt();
 721     eta = track->Eta();
 722     phi = track->Phi();
 723
 724     fhPtNoCut->Fill(pt);
 725
 726     AliAODTrack * aodTrack = dynamic_cast<AliAODTrack*>(track);
 727     AliESDtrack * esdTrack = dynamic_cast<AliESDtrack*>(track);
 728
 729     //TOF
 730     ULong_t status = track->GetStatus();
 731     Bool_t okTOF = (status & AliVTrack::kTOFout) == AliVTrack::kTOFout ;
 732     Double32_t tof = track->GetTOFsignal()*1e-3;
 733
 734     //DCA
 735     Double_t dcaCons = -999;
 736     if(aodTrack)
 737     {
 738       dcaCons = aodTrack->DCA();
 739       //vtxBC   = aodTrack->GetProdVertex()->GetBC();
 740     }
 741
 742     Double_t dca[2]   = {1e6,1e6};
 743     Double_t covar[3] = {1e6,1e6,1e6};
 744     track->PropagateToDCA(GetReader()->GetInputEvent()->GetPrimaryVertex(),bz,100.,dca,covar);
 745
 746     Float_t trackDCA = dca[0];
 747
 748     if(dcaCons == -999)
 749     {
 750       fhPtDCA[0]->Fill(pt, dca[0]);
 751       fhPtDCA[1]->Fill(pt, dca[1]);
 752     }
 753     else
 754     {
 755       trackDCA = dcaCons;
 756       fhPtDCA[2]->Fill(pt, dcaCons);
 757     }
 758
 759     if(GetReader()->AcceptDCA(pt,trackDCA)) fhPtCutDCA->Fill(pt);
 760
 761     if(fFillVertexBC0Histograms)
 762     {
 763       if(TMath::Abs(vtxBC) > 0 && vtxBC!=AliVTrack::kTOFBCNA)
 764       {
 765         fhPtVtxOutBC0->Fill(pt);
 766         fhEtaPhiVtxOutBC0->Fill(eta,phi);
 767
 768         if(dcaCons == -999)
 769         {
 770           fhPtDCAVtxOutBC0[0]->Fill(pt, dca[0]);
 771           fhPtDCAVtxOutBC0[1]->Fill(pt, dca[1]);
 772         }
 773         else
 774           fhPtDCAVtxOutBC0[2]->Fill(pt, dcaCons);
 775       }
 776       else
 777       {
 778         fhPtVtxInBC0->Fill(pt);
 779         fhEtaPhiVtxInBC0->Fill(eta,phi);
 780         if(GetReader()->AcceptDCA(pt,trackDCA)) fhPtCutDCABCOK->Fill(pt);
 781
 782         if(dcaCons == -999)
 783         {
 784           fhPtDCAVtxInBC0[0]->Fill(pt, dca[0]);
 785           fhPtDCAVtxInBC0[1]->Fill(pt, dca[1]);
 786         }
 787         else
 788           fhPtDCAVtxInBC0[2]->Fill(pt, dcaCons);
 789
 790       }
 791     }
 792
 793     if(fFillPileUpHistograms && GetReader()->IsPileUpFromSPD())
 794     {
 795       if(dcaCons == -999)
 796       {
 797         fhPtDCAPileUp[0]->Fill(pt, dca[0]);
 798         fhPtDCAPileUp[1]->Fill(pt, dca[1]);
 799       }
 800       else
 801         fhPtDCAPileUp[2]->Fill(pt, dcaCons);
 802
 803       if(fFillVertexBC0Histograms)
 804       {
 805         if(TMath::Abs(vtxBC) > 0 && vtxBC!=AliVTrack::kTOFBCNA)
 806         {
 807           if(dcaCons == -999)
 808           {
 809             fhPtDCAVtxOutBC0PileUp[0]->Fill(pt, dca[0]);
 810             fhPtDCAVtxOutBC0PileUp[1]->Fill(pt, dca[1]);
 811           }
 812           else fhPtDCAVtxOutBC0PileUp[2]->Fill(pt, dcaCons);
 813         }
 814         else
 815         {
 816           if(dcaCons == -999)
 817           {
 818             fhPtDCAVtxInBC0PileUp[0]->Fill(pt, dca[0]);
 819             fhPtDCAVtxInBC0PileUp[1]->Fill(pt, dca[1]);
 820           }
 821           else fhPtDCAVtxInBC0PileUp[2]->Fill(pt, dcaCons);
 822         }
 823       }
 824     }
 825
 826     if(!okTOF)
 827     {
 828       if(dcaCons == -999)
 829       {
 830         fhPtDCANoTOFHit[0]->Fill(pt, dca[0]);
 831         fhPtDCANoTOFHit[1]->Fill(pt, dca[1]);
 832       }
 833       else
 834         fhPtDCANoTOFHit[2]->Fill(pt, dcaCons);
 835
 836       if(fFillVertexBC0Histograms)
 837       {
 838         if(TMath::Abs(vtxBC) > 0 && vtxBC!=AliVTrack::kTOFBCNA)
 839         {
 840           if(dcaCons == -999)
 841           {
 842             fhPtDCAVtxOutBC0NoTOFHit[0]->Fill(pt, dca[0]);
 843             fhPtDCAVtxOutBC0NoTOFHit[1]->Fill(pt, dca[1]);
 844           }
 845           else
 846             fhPtDCAVtxOutBC0NoTOFHit[2]->Fill(pt, dcaCons);
 847         }
 848         else
 849         {
 850           if(dcaCons == -999)
 851           {
 852             fhPtDCAVtxInBC0NoTOFHit[0]->Fill(pt, dca[0]);
 853             fhPtDCAVtxInBC0NoTOFHit[1]->Fill(pt, dca[1]);
 854           }
 855           else
 856             fhPtDCAVtxInBC0NoTOFHit[2]->Fill(pt, dcaCons);
 857
 858         }
 859       }
 860
 861       if(fFillPileUpHistograms && GetReader()->IsPileUpFromSPD())
 862       {
 863         if(dcaCons == -999)
 864         {
 865           fhPtDCAPileUpNoTOFHit[0]->Fill(pt, dca[0]);
 866           fhPtDCAPileUpNoTOFHit[1]->Fill(pt, dca[1]);
 867         }
 868         else
 869           fhPtDCAPileUpNoTOFHit[2]->Fill(pt, dcaCons);
 870
 871         if(fFillVertexBC0Histograms)
 872         {
 873           if(TMath::Abs(vtxBC) > 0 && vtxBC!=AliVTrack::kTOFBCNA)
 874           {
 875             if(dcaCons == -999)
 876             {
 877               fhPtDCAVtxOutBC0PileUpNoTOFHit[0]->Fill(pt, dca[0]);
 878               fhPtDCAVtxOutBC0PileUpNoTOFHit[1]->Fill(pt, dca[1]);
 879             }
 880             else
 881               fhPtDCAVtxOutBC0PileUpNoTOFHit[2]->Fill(pt, dcaCons);
 882           }
 883           else
 884           {
 885             if(dcaCons == -999)
 886             {
 887               fhPtDCAVtxInBC0PileUpNoTOFHit[0]->Fill(pt, dca[0]);
 888               fhPtDCAVtxInBC0PileUpNoTOFHit[1]->Fill(pt, dca[1]);
 889             }
 890             else
 891               fhPtDCAVtxInBC0PileUpNoTOFHit[2]->Fill(pt, dcaCons);
 892
 893           }
 894         }
 895       }
 896     }
 897
 898     //printf("track pT %2.2f, DCA Cons %f, DCA1 %f, DCA2 %f, TOFBC %d, oktof %d, tof %f\n",
 899     //      pt,dcaCons,dca[0],dca[1],track->GetTOFBunchCrossing(bz),okTOF, tof);
 900
 901
 902 //    if( vtxBC == 0 && trackBC !=0 && trackBC!=AliVTrack::kTOFBCNA)
 903 //      printf("TOF Signal %e, BC %d, pt %f, dca_xy %f, dca_z %f, dca_tpc %f \n", tof,trackBC, pt,dca[0],dca[1],dcaCons);
 904
 905
 906     if(okTOF)
 907     {
 908       fhTOFSignal  ->Fill(tof);
 909       fhPtTOFSignal->Fill(pt, tof);
 910       if(GetReader()->AcceptDCA(pt,trackDCA)) fhPtTOFSignalDCACut->Fill(pt, tof);
 911
 912       if(fFillVertexBC0Histograms)
 913       {
 914         if(TMath::Abs(vtxBC) > 0 && vtxBC!=AliVTrack::kTOFBCNA)
 915           fhPtTOFSignalVtxOutBC0->Fill(pt, tof);
 916         else
 917           fhPtTOFSignalVtxInBC0->Fill(pt, tof);
 918       }
 919
 920       Int_t trackBC = 1000;
 921
 922       if(fFillTrackBCHistograms)
 923       {
 924         trackBC = track->GetTOFBunchCrossing(bz);
 925
 926         if(trackBC==0)
 927         {
 928           fhTOFSignalBCOK->Fill(tof);
 929
 930           if(dcaCons == -999)
 931           {
 932             fhPtDCATOFBC0[0]->Fill(pt, dca[0]);
 933             fhPtDCATOFBC0[1]->Fill(pt, dca[1]);
 934           }
 935           else
 936             fhPtDCATOFBC0[2]->Fill(pt, dcaCons);
 937
 938           if(fFillPileUpHistograms && GetReader()->IsPileUpFromSPD())
 939           {
 940             if(dcaCons == -999)
 941             {
 942               fhPtDCAPileUpTOFBC0[0]->Fill(pt, dca[0]);
 943               fhPtDCAPileUpTOFBC0[1]->Fill(pt, dca[1]);
 944             }
 945             else
 946               fhPtDCAPileUpTOFBC0[2]->Fill(pt, dcaCons);
 947           }
 948         }
 949         else if(trackBC!=AliVTrack::kTOFBCNA)
 950         {
 951           if(dcaCons == -999)
 952           {
 953             fhPtDCATOFBCOut[0]->Fill(pt, dca[0]);
 954             fhPtDCATOFBCOut[1]->Fill(pt, dca[1]);
 955           }
 956           else
 957             fhPtDCATOFBCOut[2]->Fill(pt, dcaCons);
 958
 959         }
 960       }
 961
 962       if(fFillPileUpHistograms)
 963       {
 964         if(GetReader()->IsPileUpFromSPD())               fhPtTOFSignalPileUp[0]->Fill(pt, tof);
 965         if(GetReader()->IsPileUpFromEMCal())             fhPtTOFSignalPileUp[1]->Fill(pt, tof);
 966         if(GetReader()->IsPileUpFromSPDOrEMCal())        fhPtTOFSignalPileUp[2]->Fill(pt, tof);
 967         if(GetReader()->IsPileUpFromSPDAndEMCal())       fhPtTOFSignalPileUp[3]->Fill(pt, tof);
 968         if(GetReader()->IsPileUpFromSPDAndNotEMCal())    fhPtTOFSignalPileUp[4]->Fill(pt, tof);
 969         if(GetReader()->IsPileUpFromEMCalAndNotSPD())    fhPtTOFSignalPileUp[5]->Fill(pt, tof);
 970         if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) fhPtTOFSignalPileUp[6]->Fill(pt, tof);
 971
 972         if(fFillTrackBCHistograms)
 973         {
 974           if      (trackBC ==0)  { fhEtaPhiTOFBC0    ->Fill(eta,phi); if(fFillPileUpHistograms && GetReader()->IsPileUpFromSPD()) fhEtaPhiTOFBC0PileUpSPD    ->Fill(eta,phi); }
 975           else if (trackBC < 0)  { fhEtaPhiTOFBCPlus ->Fill(eta,phi); if(fFillPileUpHistograms && GetReader()->IsPileUpFromSPD()) fhEtaPhiTOFBCPlusPileUpSPD ->Fill(eta,phi); }
 976           else if (trackBC > 0)  { fhEtaPhiTOFBCMinus->Fill(eta,phi); if(fFillPileUpHistograms && GetReader()->IsPileUpFromSPD()) fhEtaPhiTOFBCMinusPileUpSPD->Fill(eta,phi); }
 977         }
 978
 979         if(fFillVertexBC0Histograms)
 980         {
 981           if(TMath::Abs(vtxBC) > 0 && vtxBC!=AliVTrack::kTOFBCNA)
 982           {
 983             if(GetReader()->IsPileUpFromSPD())               fhPtTOFSignalVtxOutBC0PileUp[0]->Fill(pt, tof);
 984             if(GetReader()->IsPileUpFromEMCal())             fhPtTOFSignalVtxOutBC0PileUp[1]->Fill(pt, tof);
 985             if(GetReader()->IsPileUpFromSPDOrEMCal())        fhPtTOFSignalVtxOutBC0PileUp[2]->Fill(pt, tof);
 986             if(GetReader()->IsPileUpFromSPDAndEMCal())       fhPtTOFSignalVtxOutBC0PileUp[3]->Fill(pt, tof);
 987             if(GetReader()->IsPileUpFromSPDAndNotEMCal())    fhPtTOFSignalVtxOutBC0PileUp[4]->Fill(pt, tof);
 988             if(GetReader()->IsPileUpFromEMCalAndNotSPD())    fhPtTOFSignalVtxOutBC0PileUp[5]->Fill(pt, tof);
 989             if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) fhPtTOFSignalVtxOutBC0PileUp[6]->Fill(pt, tof);
 990           }
 991           else
 992           {
 993             if(GetReader()->IsPileUpFromSPD())               fhPtTOFSignalVtxInBC0PileUp[0]->Fill(pt, tof);
 994             if(GetReader()->IsPileUpFromEMCal())             fhPtTOFSignalVtxInBC0PileUp[1]->Fill(pt, tof);
 995             if(GetReader()->IsPileUpFromSPDOrEMCal())        fhPtTOFSignalVtxInBC0PileUp[2]->Fill(pt, tof);
 996             if(GetReader()->IsPileUpFromSPDAndEMCal())       fhPtTOFSignalVtxInBC0PileUp[3]->Fill(pt, tof);
 997             if(GetReader()->IsPileUpFromSPDAndNotEMCal())    fhPtTOFSignalVtxInBC0PileUp[4]->Fill(pt, tof);
 998             if(GetReader()->IsPileUpFromEMCalAndNotSPD())    fhPtTOFSignalVtxInBC0PileUp[5]->Fill(pt, tof);
 999             if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) fhPtTOFSignalVtxInBC0PileUp[6]->Fill(pt, tof);
1000           }
1001         }
1002       }
1003     }
1004
1005     //Fill AODParticle after some selection
1006     Double_t mom[3] = {track->Px(),track->Py(),track->Pz()};
1007
1008     Bool_t in = GetFiducialCut()->IsInFiducialCut(mom,"CTS") ;
1009
1010     if(GetDebug() > 1)
1011       printf("AliAnaChargedParticles::MakeAnalysisFillAOD() - Track pt %2.2f, eta %2.2f, phi %2.2f in fiducial cut %d\n",pt,eta,phi,in);
1012
1013     //Acceptance selection
1014     if(IsFiducialCutOn() && ! in ) continue ;
1015
1016     // Momentum selection
1017     if(pt < GetMinPt() || pt > GetMaxPt()) continue;
1018
1019     if(okTOF) fhTOFSignalPtCut->Fill(tof);
1020     else
1021     {
1022       fhPtTOFStatus0    ->Fill(pt);
1023       fhEtaPhiTOFStatus0->Fill(eta,phi);
1024     }
1025
1026     Bool_t bITSRefit    = (status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit;
1027     Bool_t bConstrained = kFALSE;
1028     if     (aodTrack) bConstrained = aodTrack->IsGlobalConstrained();
1029     else if(esdTrack) bConstrained = (!esdTrack->HasPointOnITSLayer(0) && !esdTrack->HasPointOnITSLayer(1));
1030     //printf("Track %d, pt %2.2f, eta %2.2f, phi %2.2f, SPDRefit %d, refit %d, dcaCons %2.2f\n",
1031     //       i, pt, eta, phi, bConstrained, bITSRefit, dcaCons);
1032
1033     if(bConstrained)
1034     {
1035       if(bITSRefit)
1036       {
1037         fhPtNoSPDRefit->Fill(pt);
1038         if(pt < 2)fhEtaPhiNoSPDRefitPt02->Fill(eta,phi);
1039         if(pt > 3)fhEtaPhiNoSPDRefitPt3 ->Fill(eta,phi);
1040
1041         if(dcaCons == -999)
1042         {
1043           fhPtDCANoSPDRefit[0]->Fill(pt, dca[0]);
1044           fhPtDCANoSPDRefit[1]->Fill(pt, dca[1]);
1045         }
1046         else
1047           fhPtDCANoSPDRefit[2]->Fill(pt, dcaCons);
1048
1049       }
1050       else
1051       {
1052         fhPtNoSPDNoRefit->Fill(pt);
1053         if(pt < 2)fhEtaPhiNoSPDNoRefitPt02->Fill(eta,phi);
1054         if(pt > 3)fhEtaPhiNoSPDNoRefitPt3 ->Fill(eta,phi);
1055         if(dcaCons == -999)
1056         {
1057           fhPtDCANoSPDNoRefit[0]->Fill(pt, dca[0]);
1058           fhPtDCANoSPDNoRefit[1]->Fill(pt, dca[1]);
1059         }
1060         else
1061           fhPtDCANoSPDNoRefit[2]->Fill(pt, dcaCons);
1062
1063       }
1064     }
1065     else
1066     {
1067       fhPtSPDRefit->Fill(pt);
1068       if(pt < 2)fhEtaPhiSPDRefitPt02->Fill(eta,phi);
1069       if(pt > 3)fhEtaPhiSPDRefitPt3 ->Fill(eta,phi);
1070       if(dcaCons == -999)
1071       {
1072         fhPtDCASPDRefit[0]->Fill(pt, dca[0]);
1073         fhPtDCASPDRefit[1]->Fill(pt, dca[1]);
1074       }
1075       else
1076         fhPtDCASPDRefit[2]->Fill(pt, dcaCons);
1077     }
1078
1079     // Mixed event
1080     if (GetMixedEvent())
1081     {
1082       evtIndex = GetMixedEvent()->EventIndex(track->GetID()) ;
1083     }
1084
1085     GetVertex(vert,evtIndex);
1086     if(TMath::Abs(vert[2])> GetZvertexCut()) return;
1087
1088     AliAODPWG4Particle tr = AliAODPWG4Particle(mom[0],mom[1],mom[2],0);
1089     tr.SetDetector("CTS");
1090     tr.SetLabel(track->GetLabel());
1091     tr.SetTrackLabel(track->GetID(),-1);
1092     tr.SetChargedBit(track->Charge()>0);
1093
1094     AddAODParticle(tr);
1095
1096   }//loop
1097
1098   if(GetDebug() > 0)
1099     printf("AliAnaChargedParticles::MakeAnalysisFillAOD() - Final aod branch entries %d\n", GetOutputAODBranch()->GetEntriesFast());
1100 }
1101
1102 //__________________________________________________________________
1103 void  AliAnaChargedParticles::MakeAnalysisFillHistograms()
1104 {
1105   //Do analysis and fill histograms
1106
1107   //Loop on stored AODParticles
1108   Int_t naod = GetOutputAODBranch()->GetEntriesFast();
1109
1110   fhNtracks->Fill(GetReader()->GetTrackMultiplicity()) ;
1111
1112   if(GetDebug() > 0)
1113     printf("AliAnaChargedParticles::MakeAnalysisFillHistograms() - aod branch entries %d\n", naod);
1114
1115   Float_t pt  = 0;
1116   Float_t phi = 0;
1117   Float_t eta = 0;
1118
1119   for(Int_t iaod = 0; iaod < naod ; iaod++)
1120   {
1121     AliAODPWG4Particle* track =  (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
1122
1123     pt  = track->Pt();
1124     eta = track->Eta();
1125     phi = track->Phi();
1126
1127     fhPt->Fill(pt);
1128
1129     if(track->GetChargedBit())
1130     {
1131       fhPhiPos   ->Fill(pt, phi);
1132       fhEtaPos   ->Fill(pt, eta);
1133       fhEtaPhiPos->Fill(eta,phi);
1134     }
1135     else
1136     {
1137       fhPhiNeg   ->Fill(pt, phi);
1138       fhEtaNeg   ->Fill(pt, eta);
1139       fhEtaPhiNeg->Fill(eta,phi);
1140     }
1141
1142     if(fFillPileUpHistograms)
1143     {
1144       if(GetReader()->IsPileUpFromSPD())               {fhPtPileUp[0]->Fill(pt);}
1145       if(GetReader()->IsPileUpFromEMCal())             {fhPtPileUp[1]->Fill(pt);}
1146       if(GetReader()->IsPileUpFromSPDOrEMCal())        {fhPtPileUp[2]->Fill(pt);}
1147       if(GetReader()->IsPileUpFromSPDAndEMCal())       {fhPtPileUp[3]->Fill(pt);}
1148       if(GetReader()->IsPileUpFromSPDAndNotEMCal())    {fhPtPileUp[4]->Fill(pt);}
1149       if(GetReader()->IsPileUpFromEMCalAndNotSPD())    {fhPtPileUp[5]->Fill(pt);}
1150       if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) {fhPtPileUp[6]->Fill(pt);}
1151     }
1152
1153
1154     if(IsDataMC())
1155     {
1156       //Play with the MC stack if available
1157       Int_t mompdg = -1;
1158       Int_t label  = track->GetLabel();
1159
1160       if(label >= 0)
1161         {
1162           if( GetReader()->ReadStack() && label < GetMCStack()->GetNtrack())
1163             {
1164               TParticle * mom = GetMCStack()->Particle(label);
1165               mompdg =TMath::Abs(mom->GetPdgCode());
1166             }
1167           else if(GetReader()->ReadAODMCParticles())
1168             {
1169               AliAODMCParticle * aodmom = 0;
1170               //Get the list of MC particles
1171               aodmom = (AliAODMCParticle*) (GetReader()->GetAODMCParticles())->At(label);
1172               mompdg =TMath::Abs(aodmom->GetPdgCode());
1173             }
1174         }
1175
1176       if(mompdg==211 || mompdg==2212 || mompdg==321 || mompdg==11){
1177         if(TMath::Abs(eta) < 0.8 && phi < 3.145) fhMCRecPt->Fill(pt);
1178       }
1179
1180       if(mompdg==211)
1181       {
1182         fhPtPion ->Fill(pt);
1183         fhPhiPion->Fill(pt, phi);
1184         fhEtaPion->Fill(pt, eta);
1185       }
1186       else if(mompdg==2212)
1187       {
1188         fhPtProton ->Fill(pt);
1189         fhPhiProton->Fill(pt, phi);
1190         fhEtaProton->Fill(pt, eta);
1191       }
1192       else if(mompdg==321)
1193       {
1194         fhPtKaon ->Fill(pt);
1195         fhPhiKaon->Fill(pt, phi);
1196         fhEtaKaon->Fill(pt, eta);
1197       }
1198       else if(mompdg==11)
1199       {
1200         fhPtElectron ->Fill(pt);
1201         fhPhiElectron->Fill(pt, phi);
1202         fhEtaElectron->Fill(pt, eta);
1203       }
1204       else
1205       {
1206         fhPtUnknown ->Fill(pt);
1207         fhPhiUnknown->Fill(pt, phi);
1208         fhEtaUnknown->Fill(pt, eta);
1209       }
1210     }//Work with stack also
1211
1212   }// aod branch loop
1213
1214   if(IsDataMC())
1215     {
1216       Int_t primpdg = -1;
1217       if(GetReader()->ReadStack()) {
1218         AliStack * stack = GetMCStack();
1219         if(stack) {
1220           for(Int_t i=0 ; i<stack->GetNtrack(); i++){
1221             TParticle *prim = stack->Particle(i);
1222             primpdg = TMath::Abs(prim->GetPdgCode());
1223             if(primpdg == 211 || primpdg == 2212 || primpdg == 321 || primpdg == 11){
1224               if(prim->IsPrimary() && TMath::Abs(prim->Eta()) < 0.8 && prim->Phi() < 3.145){
1225                 fhMCPt->Fill(prim->Pt());
1226                 fhMCPhi->Fill(prim->Pt(),prim->Phi());
1227                 fhMCEta->Fill(prim->Pt(),prim->Eta());
1228               }
1229             }
1230           }
1231         }
1232       }
1233       else if(GetReader()->ReadAODMCParticles()){
1234         TClonesArray * mcparticles = GetReader()->GetAODMCParticles();
1235         if(mcparticles){
1236           Int_t nprim = mcparticles->GetEntriesFast();
1237           for(Int_t i=0; i < nprim; i++){
1238             AliAODMCParticle * prim = (AliAODMCParticle *) mcparticles->At(i);
1239             primpdg = TMath::Abs(prim->GetPdgCode());
1240             if(primpdg == 211 || primpdg == 2212 || primpdg == 321 || primpdg == 11){
1241               if(prim->IsPhysicalPrimary() && TMath::Abs(prim->Eta()) < 0.8 && prim->Phi() < 3.145){
1242                 fhMCPt->Fill(prim->Pt());
1243                 fhMCPhi->Fill(prim->Pt(),prim->Phi());
1244                 fhMCEta->Fill(prim->Pt(),prim->Eta());
1245               }
1246             }
1247           }
1248         }
1249       }
1250     }
1251   }