PWGGA/GammaConv/AliAnalysisTaskGCPartToPWG4Part.cxx

   1 #include <iostream>
   2 #include "TChain.h"
   3 #include "TTree.h"
   4 #include "TH1F.h"
   5 #include "TCanvas.h"
   6 #include "TString.h"
   7
   8 #include "AliAnalysisTask.h"
   9 #include "AliAnalysisManager.h"
  10 #include "AliAnalysisTaskGCPartToPWG4Part.h"
  11
  12 #include "AliESDEvent.h"
  13 #include "AliESDCaloCluster.h"
  14 #include "AliESDInputHandler.h"
  15
  16 #include "AliAODPWG4ParticleCorrelation.h"
  17 #include "AliAODEvent.h"
  18 #include "AliAODHandler.h"
  19 #include "AliAODCaloCluster.h"
  20 #include "AliAODConversionPhoton.h"
  21 #include "AliAODJet.h"
  22
  23 #include "AliAODInputHandler.h"
  24
  25 #include "AliAODMCParticle.h"
  26
  27
  28 #include "AliAODMCHeader.h"
  29 // Gamma - jet correlation analysis task
  30 // Authors: Svein Lindal
  31
  32
  33 using namespace std;
  34
  35 ClassImp(AliAnalysisTaskGCPartToPWG4Part)
  36
  37 //________________________________________________________________________
  38 AliAnalysisTaskGCPartToPWG4Part::AliAnalysisTaskGCPartToPWG4Part()
  39 : AliAnalysisTaskSE(),
  40   fDeltaAODFileName(""),
  41   fGammaCutString("GammaConv"),
  42   fPionCutString("GammaConv"),
  43   fAODBranchName("GammaConv_gamma"),
  44   fAODPWG4Photons(NULL),
  45   fAODPWG4Pi0(NULL),
  46    fDebugLevel(0)
  47 {
  48   // Dummy Constructor
  49 }
  50
  51 //________________________________________________________________________________
  52 AliAnalysisTaskGCPartToPWG4Part::~AliAnalysisTaskGCPartToPWG4Part() {
  53
  54   // dtor
  55
  56   if(fAODPWG4Photons)
  57     delete fAODPWG4Photons;
  58   fAODPWG4Photons = NULL;
  59
  60   if(fAODPWG4Pi0)
  61     delete fAODPWG4Pi0;
  62   fAODPWG4Pi0 = NULL;
  63
  64
  65 }
  66
  67
  68
  69 //________________________________________________________________________
  70 AliAnalysisTaskGCPartToPWG4Part::AliAnalysisTaskGCPartToPWG4Part(const char *name) :
  71   AliAnalysisTaskSE(name),
  72   fDeltaAODFileName(""),
  73   fGammaCutString("GammaConv"),
  74   fPionCutString("GammaConv"),
  75   fAODBranchName("GammaConv_gamma"),
  76   fAODPWG4Photons(NULL),
  77   fAODPWG4Pi0(NULL),
  78   fDebugLevel(0)
  79 {
  80   // Constructor
  81
  82   DefineInput(0, TChain::Class());
  83   // Output slot #0 id reserved by the base class for AOD
  84
  85   // Output slot #1 writes into a TH1 container
  86   DefineOutput(1, TList::Class());
  87
  88 }
  89
  90
  91
  92 //________________________________________________________________________
  93 void AliAnalysisTaskGCPartToPWG4Part::UserCreateOutputObjects() {
  94
  95   // Create branches
  96
  97   fAODPWG4Photons = new TClonesArray("AliAODPWG4ParticleCorrelation", 0);
  98   fAODPWG4Photons->SetName("PhotonsCTS");
  99   AddAODBranch("TClonesArray", &fAODPWG4Photons);
 100
 101   fAODPWG4Pi0 = new TClonesArray("AliAODPWG4ParticleCorrelation", 0);
 102   fAODPWG4Pi0->SetName("Pi0sCTS");
 103   AddAODBranch("TClonesArray", &fAODPWG4Pi0);
 104
 105 }
 106
 107 //________________________________________________________________________
 108 void AliAnalysisTaskGCPartToPWG4Part::UserExec(Option_t *)
 109 {
 110   // Main
 111
 112   //AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(kTRUE);
 113
 114   //Clear stuff for new event
 115   CleanUp();
 116
 117   ///Get AOD event
 118   AliAODEvent * aodEvent = GetAODEvent();
 119   if(!aodEvent) {
 120     AliError("No AOD event!!");
 121     return;
 122   }
 123
 124   ProcessConvGamma(aodEvent);
 125
 126 }
 127
 128 //____________________________________________________________________________________
 129 void AliAnalysisTaskGCPartToPWG4Part::ProcessConvGamma( const AliAODEvent * aodEvent )
 130 {
 131
 132   // Main method to put conversion photons/pi0 to AliAODPWG4Particle
 133
 134   TClonesArray * convGamma = GetConversionGammas(aodEvent);
 135   if(!convGamma) {
 136     AliError(Form("No gamma branch by name %s found in file %s", fAODBranchName.Data(), fDeltaAODFileName.Data()));
 137     return;
 138   }
 139
 140   //TClonesArray * arrayMC = dynamic_cast<TClonesArray*>(aodEvent->GetList()->FindObject(AliAODMCParticle::StdBranchName()));
 141   for (Int_t iPhot = 0; iPhot < convGamma->GetEntriesFast(); iPhot++) {
 142
 143     AliAODPWG4ParticleCorrelation * photon = NULL;
 144     AliAODConversionPhoton * convParticle = dynamic_cast<AliAODConversionPhoton*>(convGamma->At(iPhot));
 145     if (convParticle && BothTracksPresent(convParticle)) {
 146       photon = AddToAOD(convParticle, fAODPWG4Photons, "ConvGamma");
 147
 148     } else {
 149       continue;
 150     }
 151
 152     if(photon && fDebugLevel > 2) {
 153       printf("Added conversion photon number %d, pt: %f \n", iPhot, photon->Pt());
 154     }
 155
 156   }
 157
 158 /*
 159   TClonesArray * pions = GetPions(aodEvent);
 160   if(!pions) {
 161     AliError(Form("No pi0 branch by name %s found in file %s", fAODBranchName.Data(), fDeltaAODFileName.Data()));
 162     return;
 163   }
 164
 165   for (Int_t iPhot = 0; iPhot < pions->GetEntriesFast(); iPhot++) {
 166     AliAODPWG4ParticleCorrelation * pion = NULL;
 167     AliAODConversionPhoton * convParticle = dynamic_cast<AliAODConversionPhoton*>(pions->At(iPhot));
 168     if (convParticle && BothGammaPresent(convParticle, convGamma)) {
 169       pion = AddPionToAOD(convParticle, fAODPWG4Pi0, "ConvGamma");
 170
 171     } else {
 172       continue;
 173     }
 174
 175     if(pion && fDebugLevel > 2) {
 176       printf("Added conversion pion number %d, pt: %f \n", iPhot, pion->Pt());
 177     }
 178   }
 179 */
 180 }
 181
 182 //__________________________________________________________________________________
 183 AliAODPWG4ParticleCorrelation * AliAnalysisTaskGCPartToPWG4Part::AddToAOD(AliAODConversionPhoton * aodO, TClonesArray * branch, TString detector) {
 184   // Add photon to AOD
 185
 186   new((*branch)[branch->GetEntriesFast()]) AliAODPWG4ParticleCorrelation(aodO->Px(), aodO->Py(), aodO->Pz(), aodO->E());
 187   AliAODPWG4ParticleCorrelation * photon = dynamic_cast<AliAODPWG4ParticleCorrelation*>(branch->Last());
 188   if(photon) {
 189     photon->SetTrackLabel(aodO->GetLabel1(), aodO->GetLabel2());
 190     photon->SetDetector(detector);
 191     return photon;
 192   } else {
 193     return NULL;
 194   }
 195
 196 }
 197
 198
 199 ///__________________________________________________________________________________
 200 AliAODPWG4ParticleCorrelation * AliAnalysisTaskGCPartToPWG4Part::AddPionToAOD(AliAODConversionPhoton * pion, TClonesArray * branch, TString detector) {
 201   // Add pions to AOD
 202
 203   new((*branch)[branch->GetEntriesFast()]) AliAODPWG4ParticleCorrelation(pion->Px(), pion->Py(), pion->Pz(), pion->E());
 204   AliAODPWG4ParticleCorrelation * pwg4Pion = dynamic_cast<AliAODPWG4ParticleCorrelation*>(branch->Last());
 205   if(pwg4Pion) {
 206     Int_t tl[4] = {-1, -1, -1, -1};
 207     //pion->GetGrandChildren(photons, tl);
 208     pwg4Pion->SetTrackLabel(tl[0], tl[1], tl[2], tl[3]);
 209     pwg4Pion->SetDetector(detector);
 210     for(Int_t i = 0; i < 4; i++) {
 211       cout << tl[i] << " ";
 212     }
 213     cout << "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!"<<endl;
 214     return pwg4Pion;
 215   } else {
 216     return NULL;
 217   }
 218
 219 }
 220
 221
 222 //_____________________________________________________________________
 223 void AliAnalysisTaskGCPartToPWG4Part::Terminate(Option_t *) {
 224   // Draw result to the screen
 225   // Called once at the end of the query
 226 }
 227
 228
 229 //_____________________________________________________________________
 230 AliAODEvent * AliAnalysisTaskGCPartToPWG4Part::GetAODEvent() {
 231   //Get the AOD event from whereever it might be
 232   AliAODEvent * aodEvent = dynamic_cast<AliAODEvent*>(InputEvent());
 233   if(!aodEvent) {
 234     aodEvent = AODEvent();
 235   }
 236
 237   return aodEvent;
 238
 239 }
 240
 241
 242 //_____________________________________________________________________
 243 TClonesArray * AliAnalysisTaskGCPartToPWG4Part::GetAODBranch(const AliAODEvent * aodEvent, TString branchName) const {
 244
 245   //Get Conversion gamma branch of AOD. First try standard AOD
 246   TClonesArray * branch = dynamic_cast<TClonesArray*>(aodEvent->FindListObject(branchName.Data()));
 247
 248
 249   //If it's there, send it back
 250   if(branch)  return branch;
 251
 252
 253   //If AOD not in standard file have to locate it in delta AOD
 254   if( !(fDeltaAODFileName.Length() > 0)  ) return NULL;
 255
 256   AliAODHandler * aodHandler = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
 257   if(aodHandler) {
 258     AliAODExtension * gExt = dynamic_cast<AliAODExtension*>(aodHandler->GetExtensions()->FindObject(fDeltaAODFileName));
 259     if(gExt) {
 260       AliAODEvent * gcEvent = gExt->GetAOD();
 261       return dynamic_cast<TClonesArray*>(gcEvent->FindListObject(branchName.Data()));
 262     }
 263   }
 264   return NULL;
 265
 266 }
 267
 268
 269 //_____________________________________________________________________
 270 TClonesArray * AliAnalysisTaskGCPartToPWG4Part::GetConversionGammas(const AliAODEvent * aodEvent) const {
 271   // Returns photon branch
 272   return GetAODBranch(aodEvent, Form("%s_gamma", fGammaCutString.Data()));
 273
 274 }
 275
 276
 277 //_____________________________________________________________________
 278 TClonesArray * AliAnalysisTaskGCPartToPWG4Part::GetPions(const AliAODEvent * aodEvent) const {
 279     // returns pions branch
 280     return GetAODBranch(aodEvent, Form("%s_Pi0", fPionCutString.Data()));
 281
 282 }
 283
 284
 285 //_________________________________________________________________________
 286 void AliAnalysisTaskGCPartToPWG4Part::CleanUp() {
 287   fAODPWG4Photons->Delete();
 288   fAODPWG4Pi0->Delete();
 289 }
 290
 291
 292 //______________________________________________________________________________________________
 293 Bool_t AliAnalysisTaskGCPartToPWG4Part::BothTracksPresent(const AliAODConversionPhoton * photon)
 294 {
 295   // Check presence of electron pair
 296   AliVTrack * track1 = NULL;
 297   AliVTrack * track2 = NULL;
 298   Int_t ntracks = InputEvent()->GetNumberOfTracks();
 299   for(Int_t i = 0; i < ntracks; i++) {
 300     AliVTrack * track = (AliVTrack*)InputEvent()->GetTrack(i);
 301     if(track) {
 302       if      (track->GetID() == photon->GetLabel1()) track1 = track;
 303       else if (track->GetID() == photon->GetLabel2()) track2 = track;
 304       if(track1 && track2) break;
 305     }
 306   }
 307
 308   if(track1 && track2) {
 309     return kTRUE;
 310   }
 311   cout << "AliAnalysisTaskGCPartToPWG4Part::BothTracksPresen() - Could not get both tracks!!! labels "  << photon->GetLabel1() << " " << photon->GetLabel2()  <<endl;
 312   return kFALSE;
 313
 314 }
 315
 316
 317 //___________________________________________________________________________________________
 318 Bool_t AliAnalysisTaskGCPartToPWG4Part::BothGammaPresent(const AliAODConversionPhoton * pion,
 319                                                          const TClonesArray * photons)
 320 {
 321   // Check presence of gamma
 322   AliAODConversionPhoton * photon1 = dynamic_cast<AliAODConversionPhoton*>(photons->At(pion->GetLabel1()));
 323   AliAODConversionPhoton * photon2 = dynamic_cast<AliAODConversionPhoton*>(photons->At(pion->GetLabel2()));
 324
 325   if(photon1 && photon2) {
 326     if( BothTracksPresent(photon1) &&  BothTracksPresent(photon2)) {
 327       return kTRUE;
 328     }
 329   } else {
 330     cout << "can't find both photons "<< endl;
 331   }
 332
 333   return kFALSE;
 334 }
 335