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
  39 ClassImp(AliAnaChargedParticles)
  40
  41 //__________________________________________________
  42   AliAnaChargedParticles::AliAnaChargedParticles() :
  43     AliAnaCaloTrackCorrBaseClass(),
  44     fPdg(0),
  45     fhNtracks(0),   fhPt(0),
  46     fhPhiNeg(0),    fhEtaNeg(0),
  47     fhPhiPos(0),    fhEtaPos(0),
  48     fhEtaPhiPos(0), fhEtaPhiNeg(0),
  49     //MC
  50     fhPtPion(0),    fhPhiPion(0),     fhEtaPion(0),
  51     fhPtProton(0),  fhPhiProton(0),   fhEtaProton(0),
  52     fhPtElectron(0),fhPhiElectron(0), fhEtaElectron(0),
  53     fhPtKaon(0),    fhPhiKaon(0),     fhEtaKaon(0),
  54     fhPtUnknown(0), fhPhiUnknown(0),  fhEtaUnknown(0)
  55 {
  56   //Default Ctor
  57
  58   //Initialize parameters
  59   InitParameters();
  60
  61 }
  62
  63 //_______________________________________________________
  64 TList *  AliAnaChargedParticles::GetCreateOutputObjects()
  65 {
  66   // Create histograms to be saved in output file and
  67   // store them in fOutputContainer
  68
  69
  70   TList * outputContainer = new TList() ;
  71   outputContainer->SetName("ExampleHistos") ;
  72
  73   Int_t nptbins  = GetHistogramRanges()->GetHistoPtBins(); Int_t nphibins = GetHistogramRanges()->GetHistoPhiBins(); Int_t netabins = GetHistogramRanges()->GetHistoEtaBins();
  74   Float_t ptmax  = GetHistogramRanges()->GetHistoPtMax();  Float_t phimax = GetHistogramRanges()->GetHistoPhiMax();  Float_t etamax = GetHistogramRanges()->GetHistoEtaMax();
  75   Float_t ptmin  = GetHistogramRanges()->GetHistoPtMin();  Float_t phimin = GetHistogramRanges()->GetHistoPhiMin();  Float_t etamin = GetHistogramRanges()->GetHistoEtaMin();
  76
  77   fhNtracks  = new TH1F ("hNtracks","# of tracks", 1000,0,1000);
  78   fhNtracks->SetXTitle("# of tracks");
  79   outputContainer->Add(fhNtracks);
  80
  81   fhPt  = new TH1F ("hPt","p_T distribution", nptbins,ptmin,ptmax);
  82   fhPt->SetXTitle("p_{T} (GeV/c)");
  83   outputContainer->Add(fhPt);
  84
  85   fhPhiNeg  = new TH2F ("hPhiNegative","#phi of negative charges distribution",
  86                         nptbins,ptmin,ptmax, nphibins,phimin,phimax);
  87   fhPhiNeg->SetYTitle("#phi (rad)");
  88   fhPhiNeg->SetXTitle("p_{T} (GeV/c)");
  89   outputContainer->Add(fhPhiNeg);
  90
  91   fhEtaNeg  = new TH2F ("hEtaNegative","#eta of negative charges distribution",
  92                         nptbins,ptmin,ptmax, netabins,etamin,etamax);
  93   fhEtaNeg->SetYTitle("#eta ");
  94   fhEtaNeg->SetXTitle("p_{T} (GeV/c)");
  95   outputContainer->Add(fhEtaNeg);
  96
  97   fhPhiPos  = new TH2F ("hPhiPositive","#phi of positive charges distribution",
  98                         nptbins,ptmin,ptmax, nphibins,phimin,phimax);
  99   fhPhiPos->SetYTitle("#phi (rad)");
 100   fhPhiPos->SetXTitle("p_{T} (GeV/c)");
 101   outputContainer->Add(fhPhiPos);
 102
 103   fhEtaPos  = new TH2F ("hEtaPositive","#eta of positive charges distribution",
 104                         nptbins,ptmin,ptmax, netabins,etamin,etamax);
 105   fhEtaPos->SetYTitle("#eta ");
 106   fhEtaPos->SetXTitle("p_{T} (GeV/c)");
 107   outputContainer->Add(fhEtaPos);
 108
 109   fhEtaPhiPos  = new TH2F ("hEtaPhiPositive","pt/eta/phi of positive charge",netabins,etamin,etamax, nphibins,phimin,phimax);
 110   fhEtaPhiPos->SetXTitle("#eta ");
 111   fhEtaPhiPos->SetYTitle("#phi (rad)");
 112   outputContainer->Add(fhEtaPhiPos);
 113
 114   fhEtaPhiNeg  = new TH2F ("hEtaPhiNegative","eta vs phi of negative charge",netabins,etamin,etamax, nphibins,phimin,phimax);
 115   fhEtaPhiNeg->SetXTitle("#eta ");
 116   fhEtaPhiNeg->SetYTitle("#phi (rad)");
 117   outputContainer->Add(fhEtaPhiNeg);
 118
 119   if(IsDataMC()){
 120
 121     fhPtPion  = new TH1F ("hPtMCPion","p_T distribution from #pi", nptbins,ptmin,ptmax);
 122     fhPtPion->SetXTitle("p_{T} (GeV/c)");
 123     outputContainer->Add(fhPtPion);
 124
 125     fhPhiPion  = new TH2F ("hPhiMCPion","#phi distribution from #pi",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
 126     fhPhiPion->SetXTitle("#phi (rad)");
 127     outputContainer->Add(fhPhiPion);
 128
 129     fhEtaPion  = new TH2F ("hEtaMCPion","#eta distribution from #pi",nptbins,ptmin,ptmax, netabins,etamin,etamax);
 130     fhEtaPion->SetXTitle("#eta ");
 131     outputContainer->Add(fhEtaPion);
 132
 133     fhPtProton  = new TH1F ("hPtMCProton","p_T distribution from proton", nptbins,ptmin,ptmax);
 134     fhPtProton->SetXTitle("p_{T} (GeV/c)");
 135     outputContainer->Add(fhPtProton);
 136
 137     fhPhiProton  = new TH2F ("hPhiMCProton","#phi distribution from proton",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
 138     fhPhiProton->SetXTitle("#phi (rad)");
 139     outputContainer->Add(fhPhiProton);
 140
 141     fhEtaProton  = new TH2F ("hEtaMCProton","#eta distribution from proton",nptbins,ptmin,ptmax, netabins,etamin,etamax);
 142     fhEtaProton->SetXTitle("#eta ");
 143     outputContainer->Add(fhEtaProton);
 144
 145     fhPtKaon  = new TH1F ("hPtMCKaon","p_T distribution from kaon", nptbins,ptmin,ptmax);
 146     fhPtKaon->SetXTitle("p_{T} (GeV/c)");
 147     outputContainer->Add(fhPtKaon);
 148
 149     fhPhiKaon  = new TH2F ("hPhiMCKaon","#phi distribution from kaon",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
 150     fhPhiKaon->SetXTitle("#phi (rad)");
 151     outputContainer->Add(fhPhiKaon);
 152
 153     fhEtaKaon  = new TH2F ("hEtaMCKaon","#eta distribution from kaon",nptbins,ptmin,ptmax, netabins,etamin,etamax);
 154     fhEtaKaon->SetXTitle("#eta ");
 155     outputContainer->Add(fhEtaKaon);
 156
 157     fhPtElectron  = new TH1F ("hPtMCElectron","p_T distribution from electron", nptbins,ptmin,ptmax);
 158     fhPtElectron->SetXTitle("p_{T} (GeV/c)");
 159     outputContainer->Add(fhPtElectron);
 160
 161     fhPhiElectron  = new TH2F ("hPhiMCElectron","#phi distribution from electron",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
 162     fhPhiElectron->SetXTitle("#phi (rad)");
 163     outputContainer->Add(fhPhiElectron);
 164
 165     fhEtaElectron  = new TH2F ("hEtaMCElectron","#eta distribution from electron",nptbins,ptmin,ptmax, netabins,etamin,etamax);
 166     fhEtaElectron->SetXTitle("#eta ");
 167     outputContainer->Add(fhEtaElectron);
 168
 169     fhPtUnknown  = new TH1F ("hPtMCUnknown","p_T distribution from unknown", nptbins,ptmin,ptmax);
 170     fhPtUnknown->SetXTitle("p_{T} (GeV/c)");
 171     outputContainer->Add(fhPtUnknown);
 172
 173     fhPhiUnknown  = new TH2F ("hPhiMCUnknown","#phi distribution from unknown",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
 174     fhPhiUnknown->SetXTitle("#phi (rad)");
 175     outputContainer->Add(fhPhiUnknown);
 176
 177     fhEtaUnknown  = new TH2F ("hEtaMCUnknown","#eta distribution from unknown",nptbins,ptmin,ptmax, netabins,etamin,etamax);
 178     fhEtaUnknown->SetXTitle("#eta ");
 179     outputContainer->Add(fhEtaUnknown);
 180
 181   }
 182
 183   return outputContainer;
 184
 185 }
 186
 187 //___________________________________________
 188 void AliAnaChargedParticles::InitParameters()
 189 {
 190   //Initialize the parameters of the analysis.
 191   SetOutputAODClassName("AliAODPWG4Particle");
 192   SetOutputAODName("PWG4Particle");
 193
 194   AddToHistogramsName("AnaCharged_");
 195
 196   fPdg = -1; //Select all tracks
 197
 198 }
 199
 200 //____________________________________________________________
 201 void AliAnaChargedParticles::Print(const Option_t * opt) const
 202 {
 203   //Print some relevant parameters set for the analysis
 204   if(! opt)
 205     return;
 206
 207   printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
 208   AliAnaCaloTrackCorrBaseClass::Print(" ");
 209
 210   printf("Min Pt = %3.2f\n", GetMinPt());
 211   printf("Max Pt = %3.2f\n", GetMaxPt());
 212   printf("Select clusters with pdg %d \n",fPdg);
 213
 214 }
 215
 216 //_________________________________
 217 void AliAnaChargedParticles::Init()
 218 {
 219   //Init
 220   //Do some checks
 221   if(!GetReader()->IsCTSSwitchedOn()){
 222     printf("AliAnaChargedParticles::Init() - STOP!: You want to use CTS tracks in analysis but not read!! \n!!Check the configuration file!!\n");
 223     abort();
 224   }
 225
 226 }
 227
 228 //_________________________________________________
 229 void  AliAnaChargedParticles::MakeAnalysisFillAOD()
 230 {
 231   //Do analysis and fill aods
 232   if(!GetCTSTracks() || GetCTSTracks()->GetEntriesFast() == 0) return ;
 233
 234   Int_t ntracks = GetCTSTracks()->GetEntriesFast();
 235   Double_t vert[3] = {0,0,0}; //vertex ;
 236
 237   //Some prints
 238   if(GetDebug() > 0)
 239     printf("AliAnaChargedParticles::MakeAnalysisFillAOD() - In CTS aod entries %d\n", ntracks);
 240
 241   //Fill AODParticle with CTS aods
 242   TVector3 p3;
 243   Int_t evtIndex = 0;
 244   for(Int_t i = 0; i < ntracks; i++){
 245
 246     AliVTrack * track =  (AliVTrack*) (GetCTSTracks()->At(i));
 247
 248     //Fill AODParticle after some selection
 249     Double_t mom[3] = {track->Px(),track->Py(),track->Pz()};
 250     p3.SetXYZ(mom[0],mom[1],mom[2]);
 251
 252     Bool_t in = GetFiducialCut()->IsInFiducialCut(mom,"CTS") ;
 253
 254     if(GetDebug() > 1)
 255       printf("AliAnaChargedParticles::MakeAnalysisFillAOD() - Track pt %2.2f, phi %2.2f, eta %2.2f in fiducial cut %d\n",p3.Pt(), p3.Phi(), p3.Eta(),in);
 256
 257     //Acceptance selection
 258     if(IsFiducialCutOn() && ! in ) continue ;
 259
 260     // Momentum selection
 261     if(p3.Pt() < GetMinPt() || p3.Pt() > GetMaxPt()) continue;
 262
 263     //Keep only particles identified with fPdg
 264     //Selection not done for the moment
 265     //Should be done here.
 266
 267     // Mixed event
 268     if (GetMixedEvent()){
 269       evtIndex = GetMixedEvent()->EventIndex(track->GetID()) ;
 270     }
 271     GetVertex(vert,evtIndex);
 272     if(TMath::Abs(vert[2])> GetZvertexCut()) return;
 273
 274     AliAODPWG4Particle tr = AliAODPWG4Particle(mom[0],mom[1],mom[2],0);
 275     tr.SetDetector("CTS");
 276     tr.SetLabel(track->GetLabel());
 277     tr.SetTrackLabel(track->GetID(),-1);
 278     tr.SetChargedBit(track->Charge()>0);
 279
 280     AddAODParticle(tr);
 281
 282   }//loop
 283
 284   if(GetDebug() > 0)
 285     printf("AliAnaChargedParticles::MakeAnalysisFillAOD() - Final aod branch entries %d\n", GetOutputAODBranch()->GetEntriesFast());
 286 }
 287
 288 //__________________________________________________________________
 289 void  AliAnaChargedParticles::MakeAnalysisFillHistograms()
 290 {
 291   //Do analysis and fill histograms
 292
 293   //Loop on stored AODParticles
 294   Int_t naod = GetOutputAODBranch()->GetEntriesFast();
 295
 296   fhNtracks->Fill(GetReader()->GetTrackMultiplicity()) ;
 297
 298   if(GetDebug() > 0)
 299     printf("AliAnaChargedParticles::MakeAnalysisFillHistograms() - aod branch entries %d\n", naod);
 300
 301   for(Int_t iaod = 0; iaod < naod ; iaod++){
 302     AliAODPWG4Particle* tr =  (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
 303
 304     fhPt->Fill(tr->Pt());
 305
 306     if(tr->GetChargedBit()){
 307       fhPhiPos   ->Fill(tr->Pt(), tr->Phi());
 308       fhEtaPos   ->Fill(tr->Pt(), tr->Eta());
 309       fhEtaPhiPos->Fill(tr->Eta(),tr->Phi());
 310     }
 311     else{
 312       fhPhiNeg   ->Fill(tr->Pt(), tr->Phi());
 313       fhEtaNeg   ->Fill(tr->Pt(), tr->Eta());
 314       fhEtaPhiNeg->Fill(tr->Eta(),tr->Phi());
 315     }
 316
 317     if(IsDataMC()){
 318       //Play with the MC stack if available
 319       Int_t mompdg = -1;
 320       Int_t label  = TMath::Abs(tr->GetLabel());
 321       if(GetReader()->ReadStack()){
 322         TParticle * mom = GetMCStack()->Particle(label);
 323         mompdg =TMath::Abs(mom->GetPdgCode());
 324       }
 325       else if(GetReader()->ReadAODMCParticles()){
 326         AliAODMCParticle * aodmom = 0;
 327         //Get the list of MC particles
 328         aodmom = (AliAODMCParticle*) (GetReader()->GetAODMCParticles(tr->GetInputFileIndex()))->At(label);
 329         mompdg =TMath::Abs(aodmom->GetPdgCode());
 330       }
 331
 332       if(mompdg==211){
 333         fhPtPion ->Fill(tr->Pt());
 334         fhPhiPion->Fill(tr->Pt(), tr->Phi());
 335         fhEtaPion->Fill(tr->Pt(), tr->Eta());
 336       }
 337       else if(mompdg==2212){
 338         fhPtProton ->Fill(tr->Pt());
 339         fhPhiProton->Fill(tr->Pt(), tr->Phi());
 340         fhEtaProton->Fill(tr->Pt(), tr->Eta());
 341       }
 342       else if(mompdg==321){
 343         fhPtKaon ->Fill(tr->Pt());
 344         fhPhiKaon->Fill(tr->Pt(), tr->Phi());
 345         fhEtaKaon->Fill(tr->Pt(), tr->Eta());
 346       }
 347       else if(mompdg==11){
 348         fhPtElectron ->Fill(tr->Pt());
 349         fhPhiElectron->Fill(tr->Pt(), tr->Phi());
 350         fhEtaElectron->Fill(tr->Pt(), tr->Eta());
 351       }
 352       else {
 353         //printf("unknown pdg %d\n",mompdg);
 354         fhPtUnknown ->Fill(tr->Pt());
 355         fhPhiUnknown->Fill(tr->Pt(), tr->Phi());
 356         fhEtaUnknown->Fill(tr->Pt(), tr->Eta());
 357       }
 358
 359     }//Work with stack also
 360
 361   }// aod branch loop
 362
 363 }