PWGGA/CaloTrackCorrelations/AliAnaParticlePartonCorrelation.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 // Class for the analysis of particle-parton correlations
  18 // Particle (for example direct gamma) must be found in a previous analysis
  19 // -- Author: Gustavo Conesa (LNF-INFN)
  20 //////////////////////////////////////////////////////////////////////////////
  21
  22 // --- ROOT system ---
  23 //#include "Riostream.h"
  24 #include "TH2F.h"
  25 #include "TParticle.h"
  26 #include "TClass.h"
  27
  28 //---- ANALYSIS system ----
  29 #include "AliAnaParticlePartonCorrelation.h"
  30 #include "AliStack.h"
  31 #include "AliAODPWG4ParticleCorrelation.h"
  32
  33 ClassImp(AliAnaParticlePartonCorrelation)
  34
  35
  36 //________________________________________________________________
  37 AliAnaParticlePartonCorrelation::AliAnaParticlePartonCorrelation() :
  38 AliAnaCaloTrackCorrBaseClass(),
  39 fhDeltaEtaNearParton(0), fhDeltaPhiNearParton(0),
  40 fhDeltaPtNearParton(0),  fhPtRatNearParton(0),
  41 fhDeltaEtaAwayParton(0), fhDeltaPhiAwayParton(0),
  42 fhDeltaPtAwayParton(0),  fhPtRatAwayParton(0)
  43 {
  44   //Default Ctor
  45
  46   //Initialize parameters
  47   InitParameters();
  48 }
  49
  50 //________________________________________________________________
  51 TList *  AliAnaParticlePartonCorrelation::GetCreateOutputObjects()
  52 {
  53   // Create histograms to be saved in output file
  54
  55   TList * outputContainer = new TList() ;
  56   outputContainer->SetName("ParticlePartonHistos") ;
  57
  58   fhDeltaPhiNearParton  = new TH2F
  59   ("DeltaPhiNearParton","#phi_{particle} - #phi_{parton} vs p_{T particle}",
  60    200,0,120,200,0,6.4);
  61   fhDeltaPhiNearParton->SetYTitle("#Delta #phi");
  62   fhDeltaPhiNearParton->SetXTitle("p_{T particle} (GeV/c)");
  63   outputContainer->Add(fhDeltaPhiNearParton);
  64
  65   fhDeltaEtaNearParton  = new TH2F
  66   ("DeltaEtaNearParton","#eta_{particle} - #eta_{parton} vs p_{T particle}",
  67    200,0,120,200,-2,2);
  68   fhDeltaEtaNearParton->SetYTitle("#Delta #eta");
  69   fhDeltaEtaNearParton->SetXTitle("p_{T particle} (GeV/c)");
  70   outputContainer->Add(fhDeltaEtaNearParton);
  71
  72   fhDeltaPtNearParton  = new TH2F
  73   ("DeltaPtNearParton","#p_{T particle} - #p_{T parton} vs p_{T particle}",
  74    200,0,120,100,-10,10);
  75   fhDeltaPtNearParton->SetYTitle("#Delta #p_{T}");
  76   fhDeltaPtNearParton->SetXTitle("p_{T particle} (GeV/c)");
  77   outputContainer->Add(fhDeltaPtNearParton);
  78
  79   fhPtRatNearParton  = new TH2F
  80   ("PtRatNearParton","#p_{T parton} / #p_{T particle} vs p_{T particle}",
  81    200,0,120,200,0,5);
  82   fhPtRatNearParton->SetYTitle("ratio");
  83   fhPtRatNearParton->SetXTitle("p_{T particle} (GeV/c)");
  84   outputContainer->Add(fhPtRatNearParton);
  85
  86   fhDeltaPhiAwayParton  = new TH2F
  87   ("DeltaPhiAwayParton","#phi_{particle} - #phi_{parton} vs p_{T particle}",
  88    200,0,120,200,0,6.4);
  89   fhDeltaPhiAwayParton->SetYTitle("#Delta #phi");
  90   fhDeltaPhiAwayParton->SetXTitle("p_{T particle} (GeV/c)");
  91   outputContainer->Add(fhDeltaPhiAwayParton);
  92
  93   fhDeltaEtaAwayParton  = new TH2F
  94   ("DeltaEtaAwayParton","#eta_{particle} - #eta_{parton} vs p_{T particle}",
  95    200,0,120,200,-2,2);
  96   fhDeltaEtaAwayParton->SetYTitle("#Delta #eta");
  97   fhDeltaEtaAwayParton->SetXTitle("p_{T particle} (GeV/c)");
  98   outputContainer->Add(fhDeltaEtaAwayParton);
  99
 100   fhDeltaPtAwayParton  = new TH2F
 101   ("DeltaPtAwayParton","#p_{T particle} - #p_{T parton} vs p_{T particle}",
 102    200,0,120,100,-10,10);
 103   fhDeltaPtAwayParton->SetYTitle("#Delta #p_{T}");
 104   fhDeltaPtAwayParton->SetXTitle("p_{T particle} (GeV/c)");
 105   outputContainer->Add(fhDeltaPtAwayParton);
 106
 107   fhPtRatAwayParton  = new TH2F
 108   ("PtRatAwayParton","#p_{T parton} / #p_{T particle} vs p_{T particle}",
 109    200,0,120,200,0,5);
 110   fhPtRatAwayParton->SetYTitle("ratio");
 111   fhPtRatAwayParton->SetXTitle("p_{T particle} (GeV/c)");
 112   outputContainer->Add(fhPtRatAwayParton);
 113
 114   return outputContainer;
 115
 116 }
 117
 118 //____________________________________________________
 119 void AliAnaParticlePartonCorrelation::InitParameters()
 120 {
 121
 122   //Initialize the parameters of the analysis.
 123   SetInputAODName("PWG4Particle");
 124   SetAODObjArrayName("Partons");
 125   AddToHistogramsName("AnaPartonCorr_");
 126
 127 }
 128
 129 //_____________________________________________________________________
 130 void AliAnaParticlePartonCorrelation::Print(const Option_t * opt) const
 131 {
 132
 133   //Print some relevant parameters set for the analysis
 134   if(! opt)
 135     return;
 136
 137   printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
 138   AliAnaCaloTrackCorrBaseClass::Print(" ");
 139
 140 }
 141
 142 //__________________________________________________________
 143 void  AliAnaParticlePartonCorrelation::MakeAnalysisFillAOD()
 144 {
 145   //Particle-Parton Correlation Analysis, create AODs
 146   //Add partons to the reference list of the trigger particle
 147   //Partons are considered those in the first eight possitions in the stack
 148   //being 0, and 1 the 2 protons, and 6 and 7 the outgoing final partons.
 149   if(!GetInputAODBranch()){
 150     printf("AliAnaParticlePartonCorrelation::MakeAnalysisFillAOD() - No input particles in AOD with name branch < %s > \n",GetInputAODName().Data());
 151     abort();
 152   }
 153
 154   if(strcmp(GetInputAODBranch()->GetClass()->GetName(), "AliAODPWG4ParticleCorrelation")){
 155     printf("AliAnaParticlePartonCorrelation::MakeAnalysisFillAOD() - Wrong type of AOD object, change AOD class name in input AOD: It should be <AliAODPWG4ParticleCorrelation> and not <%s> \n",GetInputAODBranch()->GetClass()->GetName());
 156     abort();
 157   }
 158
 159   if(GetDebug() > 1){
 160     printf("AliAnaParticlePartonCorrelation::MakeAnalysisFillAOD() - Begin fill AODs \n");
 161     printf("AliAnaParticlePartonCorrelation::MakeAnalysisFillAOD() - In particle branch aod entries %d\n", GetInputAODBranch()->GetEntriesFast());
 162   }
 163
 164   //Loop on stored AOD particles
 165   Int_t naod = GetInputAODBranch()->GetEntriesFast();
 166   for(Int_t iaod = 0; iaod < naod ; iaod++){
 167     AliAODPWG4ParticleCorrelation* particle =  (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
 168
 169     AliStack * stack =  GetMCStack() ;
 170     if(!stack){
 171       printf("AliAnaParticlePartonCorrelation::MakeAnalysisFillAOD() - No Stack available, STOP\n");
 172       abort();
 173     }
 174     if(stack->GetNtrack() < 8) {
 175       printf("AliAnaParticlePartonCorrelation::MakeAnalysisFillAOD() *** small number of particles, not a PYTHIA simulation? ***:  n tracks %d \n", stack->GetNprimary());
 176       continue ;
 177     }
 178
 179     //Fill AOD reference only with partons
 180
 181     //Array with reference to partons, initialize
 182     TObjArray * objarray  = NULL;
 183     Int_t nrefs = 0;
 184
 185     TParticle * parton    = NULL ;
 186     for(Int_t ipr = 0;ipr < 8; ipr ++ ){
 187       parton = stack->Particle(ipr) ;
 188       nrefs++;
 189       if(nrefs==1){
 190         objarray = new TObjArray(0);
 191         objarray->SetName(GetAODObjArrayName());
 192         objarray->SetOwner(kFALSE);
 193       }
 194       objarray->Add(parton);
 195     }//parton loop
 196
 197     if(objarray->GetEntriesFast() > 0) particle->AddObjArray(objarray);
 198
 199   }//Aod branch loop
 200
 201   if(GetDebug() > 1) printf("AliAnaParticlePartonCorrelation::MakeAnalysisFillAOD() - End fill AODs \n");
 202 }
 203
 204 //_________________________________________________________________
 205 void  AliAnaParticlePartonCorrelation::MakeAnalysisFillHistograms()
 206 {
 207   //Particle-Parton Correlation Analysis, fill histograms
 208   if(!GetInputAODBranch()){
 209     printf("AliAnaParticlePartonCorrelation::MakeAnalysisFillHistograms() - No input particles in AOD with name branch < %s > \n",GetInputAODName().Data());
 210     abort();
 211   }
 212   if(GetDebug() > 1){
 213     printf("AliAnaParticlePartonCorrelation::MakeAnalysisFillHistograms() - Begin parton correlation analysis, fill histograms \n");
 214     printf("AliAnaParticlePartonCorrelation::MakeAnalysisFillHistograms() - In particle branch aod entries %d\n", GetInputAODBranch()->GetEntriesFast());
 215   }
 216
 217   AliStack * stack =  GetMCStack() ;
 218   if(!stack) {
 219     printf("AliAnaParticlePartonCorrelation::MakeAnalysisFillHistograms() - No Stack available, STOP\n");
 220     abort();
 221   }
 222
 223   //Loop on stored AOD particles
 224   Int_t naod = GetInputAODBranch()->GetEntriesFast();
 225   TParticle *  mom = NULL ;
 226
 227   for(Int_t iaod = 0; iaod < naod ; iaod++){
 228     AliAODPWG4ParticleCorrelation* particle =  (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
 229
 230     Float_t ptTrigg  = particle->Pt();
 231     Float_t phiTrigg = particle->Phi();
 232     Float_t etaTrigg = particle->Eta();
 233     Int_t imom = particle->GetLabel();
 234     Int_t iparent  = 2000;
 235     Int_t iawayparent = -1;
 236
 237     TObjArray * objarray = particle->GetObjArray(GetAODObjArrayName());
 238     if(!(objarray) || (objarray->GetEntriesFast() < 7) ) {
 239       printf("AliAnaParticlePartonCorrelation::MakeAnalysisFillHistograms() - Reference list with partons not filled, STOP analysis\n");
 240       abort();
 241     }
 242
 243     //Check and get indeces of mother and parton
 244     if(imom < 8 ) iparent = imom ;   //mother is already a parton
 245     else if (imom <  stack->GetNtrack()) {
 246       mom =  stack->Particle(imom);
 247       if(mom){
 248         iparent=mom->GetFirstMother();
 249         //cout<<" iparent "<<iparent<<endl;
 250         while(iparent > 7 ){
 251           mom = stack->Particle(iparent);
 252           if (mom) {
 253             imom = iparent ; //Mother label is of the inmediate parton daughter
 254             iparent = mom->GetFirstMother();
 255           }
 256           else iparent = -1;
 257           //cout<<" while iparent "<<iparent<<endl;
 258         }
 259       }
 260     }
 261
 262     if(GetDebug() > 1) printf("AliAnaParticlePartonCorrelation::MakeAnalysisFillHistograms() - N reference partons %d; labels:  mother %d, parent %d \n", objarray->GetEntriesFast(), imom, iparent);
 263
 264
 265     if(iparent < 0 || iparent > 8) {
 266       if(GetDebug() > 0 ) printf("AliAnaParticlePartonCorrelation::MakeAnalysisFillHistograms() - Failed to find appropriate parton, index %d", iparent);
 267       continue ;
 268     }
 269
 270     //Near parton is the parton that fragmented and created the mother
 271     TParticle * nearParton = (TParticle*) objarray->At(iparent);
 272     Float_t  ptNearParton    = nearParton->Pt();
 273     Float_t  phiNearParton   = nearParton->Phi() ;
 274     Float_t  etaNearParton   = nearParton->Eta() ;
 275
 276     fhDeltaEtaNearParton->Fill(ptTrigg,etaTrigg-etaNearParton);
 277     fhDeltaPhiNearParton->Fill(ptTrigg,phiTrigg-phiNearParton);
 278     fhDeltaPtNearParton->Fill(ptTrigg,ptTrigg-ptNearParton);
 279     fhPtRatNearParton->Fill(ptTrigg,ptNearParton/ptTrigg);
 280
 281     if(iparent == 7) iawayparent =6;
 282     else if(iparent == 6) iawayparent =7;
 283     else{
 284       printf("AliAnaParticlePartonCorrelation::MakeAnalysisFillHistograms() - Parent parton is not final state, skip \n");
 285       continue;
 286     }
 287
 288     //Away parton is the other final parton.
 289     TParticle * awayParton = (TParticle*) objarray->At(iawayparent);
 290     Float_t  ptAwayParton    = awayParton->Pt();
 291     Float_t  phiAwayParton   = awayParton->Phi() ;
 292     Float_t  etaAwayParton   = awayParton->Eta() ;
 293     fhDeltaEtaAwayParton->Fill(ptTrigg,etaTrigg-etaAwayParton);
 294     fhDeltaPhiAwayParton->Fill(ptTrigg,phiTrigg-phiAwayParton);
 295     fhDeltaPtAwayParton->Fill(ptTrigg,ptTrigg-ptAwayParton);
 296     fhPtRatAwayParton->Fill(ptTrigg,ptAwayParton/ptTrigg);
 297
 298   }
 299
 300   if(GetDebug() > 1) printf("AliAnaParticlePartonCorrelation::MakeAnalysisFillHistograms() - End fill histograms \n");
 301
 302 }