PWG4/AliAnaGammaHadron.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 /* History of cvs commits:
  18  *
  19  * $Log$
  20  * Revision 1.7  2007/10/29 13:48:42  gustavo
  21  * Corrected coding violations
  22  *
  23  * Revision 1.5  2007/08/17 12:40:04  schutz
  24  * New analysis classes by Gustavo Conesa
  25  *
  26  * Revision 1.4.4.2  2007/07/26 10:32:09  schutz
  27  * new analysis classes in the the new analysis framework
  28  *
  29  *
  30  */
  31
  32 //_________________________________________________________________________
  33 // Class for the analysis of gamma - hadron correlations
  34 //*-- Author: Gustavo Conesa (LNF-INFN)
  35 //////////////////////////////////////////////////////////////////////////////
  36
  37
  38 // --- ROOT system ---
  39 #include "Riostream.h"
  40
  41 //---- AliRoot system ----
  42 #include "AliLog.h"
  43 #include "AliNeutralMesonSelection.h"
  44 #include "AliAnaGammaHadron.h"
  45
  46 ClassImp(AliAnaGammaHadron)
  47
  48
  49 //____________________________________________________________________________
  50   AliAnaGammaHadron::AliAnaGammaHadron() :
  51     AliAnaGammaCorrelation(),
  52     fhPhiCharged(0), fhPhiNeutral(0), fhEtaCharged(0), fhEtaNeutral(0),
  53     fhDeltaPhiGammaCharged(0),  fhDeltaPhiGammaNeutral(0),
  54     fhDeltaEtaGammaCharged(0), fhDeltaEtaGammaNeutral(0),
  55     fhDeltaPhiChargedPt(0),
  56     fhCorrelationGammaNeutral(0), fhCorrelationGammaCharged(0)
  57 {
  58   //Default Ctor
  59
  60   SetCorrelationType(kHadron);
  61   //Initialize parameters
  62   InitParameters();
  63 }
  64
  65 //____________________________________________________________________________
  66 AliAnaGammaHadron::AliAnaGammaHadron(const AliAnaGammaHadron & g) :
  67   AliAnaGammaCorrelation(g),
  68   fhPhiCharged(g.fhPhiCharged), fhPhiNeutral(g.fhPhiNeutral),
  69   fhEtaCharged(g.fhEtaCharged), fhEtaNeutral(g.fhEtaNeutral),
  70   fhDeltaPhiGammaCharged(g.fhDeltaPhiGammaCharged),
  71   fhDeltaPhiGammaNeutral(g.fhDeltaPhiGammaNeutral),
  72   fhDeltaEtaGammaCharged(g.fhDeltaEtaGammaCharged),
  73   fhDeltaEtaGammaNeutral(g.fhDeltaEtaGammaNeutral),
  74   fhDeltaPhiChargedPt(g.fhDeltaPhiChargedPt),
  75   fhCorrelationGammaNeutral(g.fhCorrelationGammaNeutral),
  76   fhCorrelationGammaCharged(g.fhCorrelationGammaCharged)
  77 {
  78   // cpy ctor
  79
  80 }
  81
  82 //_________________________________________________________________________
  83 AliAnaGammaHadron & AliAnaGammaHadron::operator = (const AliAnaGammaHadron & source)
  84 {
  85   // assignment operator
  86
  87   if(this == &source)return *this;
  88   ((AliAnaGammaCorrelation *)this)->operator=(source);
  89
  90   fhPhiCharged = source.fhPhiCharged ; fhPhiNeutral = source.fhPhiNeutral ;
  91   fhEtaCharged = source.fhEtaCharged ; fhEtaNeutral = source.fhEtaNeutral ;
  92   fhDeltaPhiGammaCharged = source.fhDeltaPhiGammaCharged ;
  93   fhDeltaPhiGammaNeutral = source.fhDeltaPhiGammaNeutral ;
  94   fhDeltaEtaGammaCharged = source.fhDeltaEtaGammaCharged ;
  95   fhDeltaEtaGammaNeutral = source.fhDeltaEtaGammaNeutral ;
  96   fhDeltaPhiChargedPt = source.fhDeltaPhiChargedPt ;
  97
  98   fhCorrelationGammaNeutral = source.fhCorrelationGammaNeutral ;
  99   fhCorrelationGammaCharged = source.fhCorrelationGammaCharged ;
 100
 101   return *this;
 102
 103 }
 104
 105 //____________________________________________________________________________
 106 AliAnaGammaHadron::~AliAnaGammaHadron()
 107 {
 108   // Remove all pointers except analysis output pointers.
 109 }
 110
 111
 112 //________________________________________________________________________
 113 TList *  AliAnaGammaHadron::GetCreateOutputObjects()
 114 {
 115
 116   // Create histograms to be saved in output file and
 117   // store them in fOutputContainer
 118   TList * outputContainer = new TList() ;
 119   outputContainer->SetName("GammaCorrelationHistos") ;
 120
 121   fhPhiCharged  = new TH2F
 122     ("PhiCharged","#phi_{#pi^{#pm}}  vs p_{T #gamma}",
 123      120,0,120,120,0,7);
 124   fhPhiCharged->SetYTitle("#phi_{#pi^{#pm}} (rad)");
 125   fhPhiCharged->SetXTitle("p_{T #gamma} (GeV/c)");
 126
 127   fhEtaCharged  = new TH2F
 128     ("EtaCharged","#eta_{#pi^{#pm}}  vs p_{T #gamma}",
 129      120,0,120,120,-1,1);
 130   fhEtaCharged->SetYTitle("#eta_{#pi^{#pm}} (rad)");
 131   fhEtaCharged->SetXTitle("p_{T #gamma} (GeV/c)");
 132
 133   fhDeltaPhiGammaCharged  = new TH2F
 134     ("DeltaPhiGammaCharged","#phi_{#gamma} - #phi_{charged #pi} vs p_{T #gamma}",
 135      200,0,120,200,0,6.4);
 136   fhDeltaPhiGammaCharged->SetYTitle("#Delta #phi");
 137   fhDeltaPhiGammaCharged->SetXTitle("p_{T #gamma} (GeV/c)");
 138
 139   fhDeltaPhiChargedPt  = new TH2F
 140     ("DeltaPhiChargedPt","#phi_{#gamma} - #phi_{charged #pi} vs p_{T #pi}",
 141      200,0,120,200,0,6.4);
 142   fhDeltaPhiChargedPt->SetYTitle("#Delta #phi");
 143   fhDeltaPhiChargedPt->SetXTitle("p_{T #pi} (GeV/c)");
 144
 145   fhDeltaEtaGammaCharged  = new TH2F
 146     ("DeltaEtaGammaCharged","#eta_{#gamma} - #eta_{#pi^{#pm}} vs p_{T #gamma}",
 147      200,0,120,200,-2,2);
 148   fhDeltaEtaGammaCharged->SetYTitle("#Delta #eta");
 149   fhDeltaEtaGammaCharged->SetXTitle("p_{T #gamma} (GeV/c)");
 150
 151   fhCorrelationGammaCharged  =
 152     new TH2F("CorrelationGammaCharged","z_{#gamma #pi} = p_{T #pi^{#pm}} / p_{T #gamma}",
 153              240,0.,120.,1000,0.,1.2);
 154   fhCorrelationGammaCharged->SetYTitle("z_{#gamma #pi}");
 155   fhCorrelationGammaCharged->SetXTitle("p_{T #gamma}");
 156
 157   outputContainer->Add(fhPhiCharged) ;
 158   outputContainer->Add(fhEtaCharged) ;
 159   outputContainer->Add(fhDeltaPhiGammaCharged) ;
 160   outputContainer->Add(fhDeltaEtaGammaCharged) ;
 161   outputContainer->Add(fhCorrelationGammaCharged) ;
 162   outputContainer->Add(fhDeltaPhiChargedPt) ;
 163
 164   if(!AreJetsOnlyInCTS()){
 165     //---- kHadron and kJetLeadCone ----
 166     fhPhiNeutral  = new TH2F
 167       ("PhiNeutral","#phi_{#pi^{0}}  vs p_{T #gamma}",
 168        120,0,120,120,0,7);
 169     fhPhiNeutral->SetYTitle("#phi_{#pi^{0}} (rad)");
 170     fhPhiNeutral->SetXTitle("p_{T #gamma} (GeV/c)");
 171
 172     fhEtaNeutral  = new TH2F
 173       ("EtaNeutral","#eta_{#pi^{0}}  vs p_{T #gamma}",
 174        120,0,120,120,-1,1);
 175     fhEtaNeutral->SetYTitle("#eta_{#pi^{0}} (rad)");
 176     fhEtaNeutral->SetXTitle("p_{T #gamma} (GeV/c)");
 177
 178     fhDeltaPhiGammaNeutral  = new TH2F
 179       ("DeltaPhiGammaNeutral","#phi_{#gamma} - #phi_{#pi^{0}} vs p_{T #gamma}",
 180        200,0,120,200,0,6.4);
 181     fhDeltaPhiGammaNeutral->SetYTitle("#Delta #phi");
 182     fhDeltaPhiGammaNeutral->SetXTitle("p_{T #gamma} (GeV/c)");
 183
 184     fhDeltaEtaGammaNeutral  = new TH2F
 185       ("DeltaEtaGammaNeutral","#eta_{#gamma} - #eta_{#pi^{#pm}} vs p_{T #gamma}",
 186        200,0,120,200,-2,2);
 187     fhDeltaEtaGammaNeutral->SetYTitle("#Delta #eta");
 188     fhDeltaEtaGammaNeutral->SetXTitle("p_{T #gamma} (GeV/c)");
 189
 190     fhCorrelationGammaNeutral  =
 191       new TH2F("CorrelationGammaNeutral","z_{#gamma #pi} = p_{T #pi^{0}} / p_{T #gamma}",
 192                240,0.,120.,1000,0.,1.2);
 193     fhCorrelationGammaNeutral->SetYTitle("z_{#gamma #pi}");
 194     fhCorrelationGammaNeutral->SetXTitle("p_{T #gamma}");
 195
 196     outputContainer->Add(fhPhiNeutral) ;
 197     outputContainer->Add(fhEtaNeutral) ;
 198     outputContainer->Add(fhDeltaPhiGammaNeutral) ;
 199     outputContainer->Add(fhDeltaEtaGammaNeutral) ;
 200     outputContainer->Add(fhCorrelationGammaNeutral) ;
 201   }
 202
 203   SetOutputContainer(outputContainer);
 204
 205   return outputContainer;
 206 }
 207
 208  //____________________________________________________________________________
 209 void AliAnaGammaHadron::InitParameters()
 210 {
 211
 212   //Initialize the parameters of the analysis.
 213
 214   SetMinPtHadron(2.)   ;
 215   SetDeltaPhiCutRange(1.5,4.5);
 216   SetJetsOnlyInCTS(kFALSE) ;
 217
 218 }
 219
 220 //__________________________________________________________________
 221 void AliAnaGammaHadron::Print(const Option_t * opt) const
 222 {
 223
 224   //Print some relevant parameters set for the analysis
 225   if(! opt)
 226     return;
 227
 228   Info("Print", "%s %s", GetName(), GetTitle() ) ;
 229   printf("Correlation analysis           =     %d\n", kHadron) ;
 230   printf("pT Hadron       >    %f\n", GetMinPtHadron()) ;
 231   printf("Phi gamma-Hadron      <     %f\n", GetDeltaPhiMaxCut()) ;
 232   printf("Phi gamma-Hadron      >     %f\n", GetDeltaPhiMinCut()) ;
 233
 234
 235 }
 236
 237 //____________________________________________________________________________
 238 void  AliAnaGammaHadron::MakeGammaCorrelation(TParticle * pGamma, TClonesArray * plCTS, TClonesArray * plCalo)
 239 {
 240   //Gamma Hadron Correlation Analysis
 241   AliDebug(2, "Make gamma-hadron correlation");
 242
 243   MakeGammaChargedCorrelation(pGamma, plCTS);
 244   if(!AreJetsOnlyInCTS())
 245   MakeGammaNeutralCorrelation(pGamma, plCalo);
 246
 247
 248 }
 249
 250 //____________________________________________________________________________
 251 void  AliAnaGammaHadron::MakeGammaChargedCorrelation(TParticle * pGamma, TClonesArray * pl)
 252 {
 253   //Gamma Charged Hadron Correlation Analysis
 254   AliDebug(2,"Make gamma-charged hadron correlation");
 255
 256   Double_t ptg  = pGamma->Pt();
 257   Double_t phig = pGamma->Phi();
 258   Double_t pt    = -100.;
 259   Double_t rat   = -100.;
 260   Double_t phi   = -100. ;
 261
 262   for(Int_t ipr = 0;ipr < pl->GetEntries() ; ipr ++ ){
 263
 264     TParticle * particle = dynamic_cast<TParticle *>(pl->At(ipr)) ;
 265
 266     pt    = particle->Pt();
 267     rat   = pt/ptg ;
 268     phi   = particle->Phi() ;
 269
 270     AliDebug(3,Form("pt %f, phi %f, phi gamma %f. Cuts:  delta phi min %f,  max%f, pT min %f",pt,phi,phig,GetDeltaPhiMinCut(),GetDeltaPhiMaxCut(),GetMinPtHadron()));
 271
 272     fhEtaCharged->Fill(ptg,particle->Eta());
 273     fhPhiCharged->Fill(ptg,phi);
 274     fhDeltaEtaGammaCharged->Fill(ptg,pGamma->Eta()-particle->Eta());
 275     fhDeltaPhiGammaCharged->Fill(ptg,phig-phi);
 276     fhDeltaPhiChargedPt->Fill(pt,phig-phi);
 277
 278     //Selection within angular and energy limits
 279     if(((phig-phi)> GetDeltaPhiMinCut()) && ((phig-phi)<GetDeltaPhiMaxCut()) && pt > GetMinPtHadron()){
 280       AliDebug(2,Form("Selected: pt %f, phi %f",pt,phi));
 281       fhCorrelationGammaCharged->Fill(ptg,rat);
 282     }
 283   }//particle loop
 284 }
 285
 286 //____________________________________________________________________________
 287 void  AliAnaGammaHadron::MakeGammaNeutralCorrelation(TParticle * pGamma, TClonesArray * pl)
 288 {
 289   //Gamma Neutral Hadron Correlation Analysis
 290   AliDebug(2,"Make gamma-neutral hadron correlation");
 291
 292   Double_t pt = -100.;
 293   Double_t rat = -100.;
 294   Double_t phi = -100. ;
 295   Double_t ptg  = pGamma->Pt();
 296   Double_t phig = pGamma->Phi();
 297
 298   TIter next(pl);
 299   TParticle * particlei = 0;
 300   TParticle * particlej = 0;
 301   TLorentzVector gammai;
 302   TLorentzVector gammaj;
 303
 304   Int_t iPrimary = -1;
 305   Int_t ksPdg = 0;
 306   Int_t jPrimary=-1;
 307
 308   while ( (particlei = (TParticle*)next()) ) {
 309     iPrimary++;
 310     ksPdg = particlei->GetPdgCode();
 311     AliDebug(2, Form("neutral particles opposite to gamma: pt %f, pdg %d", particlei->Pt(),ksPdg));
 312     //2 gamma overlapped, found with PID
 313     if(ksPdg == 111){
 314       pt  = particlei->Pt();
 315       rat = pt/ptg ;
 316       phi = particlei->Phi() ;
 317       fhEtaNeutral->Fill(ptg,particlei->Eta());
 318       fhPhiNeutral->Fill(ptg,phi);
 319       fhDeltaEtaGammaNeutral->Fill(ptg,pGamma->Eta()-particlei->Eta());
 320       fhDeltaPhiGammaNeutral->Fill(ptg,phig-phi);
 321
 322       //Selection within angular and energy limits
 323       if( (phig-phi)>GetDeltaPhiMinCut() && (phig-phi)<GetDeltaPhiMaxCut() && pt > GetMinPtHadron()){
 324         fhCorrelationGammaNeutral ->Fill(ptg,rat);
 325         AliDebug(2,Form("Selected pi0: pt %f, phi %f",pt,phi));
 326       }// cuts
 327     }// pdg = 111
 328
 329     //Make invariant mass analysis
 330     else if(ksPdg == 22){//  gamma i
 331
 332       //Search the photon companion in case it comes from  a Pi0 decay
 333       //Apply several cuts to select the good pair;
 334       particlei->Momentum(gammai);
 335        jPrimary=-1;
 336       TIter next2(pl);
 337       while ( (particlej = (TParticle*)next2()) ) {
 338         jPrimary++;
 339         if(jPrimary>iPrimary){
 340           ksPdg = particlej->GetPdgCode();
 341           particlej->Momentum(gammaj);
 342           if(ksPdg == 22 ){
 343
 344             phi = (gammai+gammaj).Phi();
 345             if(phi < 0)
 346               phi+=TMath::TwoPi();
 347             rat          = (gammai+gammaj).Pt()/ptg ;
 348
 349             //Fill histograms
 350             fhEtaNeutral->Fill(ptg,(gammai+gammaj).Eta());
 351             fhPhiNeutral->Fill(ptg,phi);
 352             fhDeltaEtaGammaNeutral->Fill(ptg,pGamma->Eta()-(gammai+gammaj).Eta());
 353             fhDeltaPhiGammaNeutral->Fill(ptg,phig-phi);
 354
 355             //Select good pair (good phit, pt cuts, aperture and invariant mass)
 356             if(GetNeutralMesonSelection()->SelectPair(pGamma, gammai, gammaj)){
 357               AliDebug(2,Form("Selected gamma pair: pt %f, phi %f",(gammai+gammaj).Pt(),phi));
 358               //correlation histogram
 359               fhCorrelationGammaNeutral ->Fill(ptg,rat);
 360             }//pair selection
 361           }//if pair of gammas
 362         }//jPrimary>iPrimary
 363       }//while
 364     }// if pdg = 22
 365   }//while
 366
 367 }