PWG4/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 "AliGammaConversionAODObject.h"
  21 #include "AliAODConversionParticle.h"
  22 #include "AliAODJet.h"
  23
  24 #include "AliAODInputHandler.h"
  25
  26 #include "AliAODMCParticle.h"
  27
  28 #include "AliMCAnalysisUtils.h"
  29
  30
  31 #include "AliAODMCHeader.h"
  32 // Gamma - jet correlation analysis task
  33 // Authors: Svein Lindal
  34
  35
  36 using namespace std;
  37
  38 ClassImp(AliAnalysisTaskGCPartToPWG4Part)
  39
  40 //________________________________________________________________________
  41 AliAnalysisTaskGCPartToPWG4Part::AliAnalysisTaskGCPartToPWG4Part()
  42 : AliAnalysisTaskSE(),
  43   fDeltaAODFileName(""),
  44   fAODBranchName("GammaConv_gamma"),
  45   fAODPWG4Particles(NULL),
  46   fAnaUtils(NULL)
  47 {
  48   // Dummy Constructor
  49 }
  50
  51 //________________________________________________________________________________
  52 AliAnalysisTaskGCPartToPWG4Part::~AliAnalysisTaskGCPartToPWG4Part() {
  53
  54   if(fAODPWG4Particles)
  55     fAODPWG4Particles = NULL;
  56   delete fAODPWG4Particles;
  57
  58   if(fAnaUtils)
  59     delete fAnaUtils;
  60   fAnaUtils = NULL;
  61
  62 }
  63
  64
  65
  66 //________________________________________________________________________
  67 AliAnalysisTaskGCPartToPWG4Part::AliAnalysisTaskGCPartToPWG4Part(const char *name) :
  68   AliAnalysisTaskSE(name),
  69   fDeltaAODFileName(""),
  70   fAODBranchName("GammaConv_gamma"),
  71   fAODPWG4Particles(NULL),
  72   fAnaUtils(NULL)
  73 {
  74   // Constructor
  75   DefineInput(0, TChain::Class());
  76   // Output slot #0 id reserved by the base class for AOD
  77
  78   // Output slot #1 writes into a TH1 container
  79   DefineOutput(1, TList::Class());
  80
  81   fAnaUtils = new AliMCAnalysisUtils();
  82 }
  83
  84
  85
  86 //________________________________________________________________________
  87 void AliAnalysisTaskGCPartToPWG4Part::UserCreateOutputObjects() {
  88   fAODPWG4Particles = new TClonesArray("AliAODPWG4ParticleCorrelation", 0);
  89   fAODPWG4Particles->SetName("ConversionGamma");
  90   AddAODBranch("TClonesArray", &fAODPWG4Particles);
  91
  92 }
  93
  94 //________________________________________________________________________
  95 void AliAnalysisTaskGCPartToPWG4Part::UserExec(Option_t *)
  96 {
  97
  98   //AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(kTRUE);
  99
 100   //Clear stuff for new event
 101   CleanUp();
 102
 103   ///Get AOD event
 104   AliAODEvent * aodEvent = GetAODEvent();
 105   if(!aodEvent) {
 106     AliError("No AOD event!!");
 107     return;
 108   }
 109
 110
 111   ProcessConvGamma(aodEvent);
 112
 113
 114   //PostData(1, fOutputList);
 115
 116 }
 117
 118
 119
 120
 121 //___________________________________________________________________________________________
 122 void AliAnalysisTaskGCPartToPWG4Part::ProcessConvGamma( const AliAODEvent * const aodEvent ) {
 123
 124   TClonesArray * tracks = aodEvent->GetTracks();
 125   if(!tracks) {
 126     cout << "No tracks!!!"<<endl;
 127     return;
 128   }
 129
 130
 131   TClonesArray * convGamma = GetConversionGammas(aodEvent);
 132   if(!convGamma) {
 133     AliError(Form("No branch by name %s found in file %s", fAODBranchName.Data(), fDeltaAODFileName.Data()));
 134     return;
 135   }
 136
 137
 138   //TClonesArray * arrayMC = dynamic_cast<TClonesArray*>(aodEvent->GetList()->FindObject(AliAODMCParticle::StdBranchName()));
 139   for (Int_t iPhot = 0; iPhot < convGamma->GetEntriesFast(); iPhot++) {
 140
 141     AliGammaConversionAODObject * aodO = dynamic_cast<AliGammaConversionAODObject*>(convGamma->At(iPhot));
 142     if(aodO) {
 143       AddToAOD(aodO, fAODPWG4Particles, "ConvGamma");
 144     } else {
 145
 146       AliAODConversionParticle * convParticle = dynamic_cast<AliAODConversionParticle*>(convGamma->At(iPhot));
 147       if (convParticle) {
 148         AddToAOD(convParticle, fAODPWG4Particles, "ConvGamma");
 149       } else {
 150         AliError(Form("ERROR: Could not receive ga %d\n", iPhot));
 151         continue;
 152       }
 153     }
 154   }
 155
 156 }
 157
 158
 159 //////_________________________________________________________________________________________
 160 Int_t AliAnalysisTaskGCPartToPWG4Part::CheckTag(AliAODPWG4ParticleCorrelation * particle, TClonesArray * tracks, TClonesArray * arrayMC, AliAODMCHeader * mcHeader) {
 161
 162
 163   // for (int imc = 0; imc < arrayMC->GetEntriesFast(); imc++) {
 164   // //for (int imc = 0; imc < 20; imc++) {
 165
 166
 167
 168   //   AliAODMCParticle * mParticle = dynamic_cast<AliAODMCParticle*>(arrayMC->At(imc));
 169   //   //    cout << mParticle->GetPdgCode() << " " <<mParticle->GetStatus() << endl;
 170
 171
 172   //   //if( mParticle->GetPdgCode() == 34) cout<< "BALLE"<<endl;// && mParticle->GetMother() < 10) {
 173
 174   //   if( mParticle->GetPdgCode() == 22 ) { //&& mParticle->GetMother() < 10) {
 175   //     //AliAODMCParticle * mother = dynamic_cast<AliAODMCParticle*>(arrayMC->At(mParticle->GetMother()));
 176   //     //if (TMath::Abs(mother->GetPdgCode()) < 22 && TMath::Abs(mother->GetPdgCode()) != 11) cout  <<"KKKKKKKKK"<< mother->GetPdgCode() <<endl;
 177
 178   //     for (int id = 0; id < mParticle->GetNDaughters(); id++) {
 179   //    int idaughter = mParticle->GetDaughter(id);
 180   //    if(idaughter>0 && idaughter < arrayMC->GetEntriesFast()) {
 181   //      AliAODMCParticle * dmParticle = dynamic_cast<AliAODMCParticle*>(arrayMC->At(mParticle->GetDaughter(id)));
 182   //      cout  << dmParticle->GetPdgCode() << endl;
 183   //    }
 184   //     }
 185   //   }
 186   // }
 187
 188
 189
 190
 191
 192
 193
 194
 195
 196
 197
 198
 199   Int_t tag = 0;
 200
 201   Int_t l1 = particle->GetTrackLabel(0);
 202   Int_t l2 = particle->GetTrackLabel(1);
 203
 204   AliAODTrack * track1 = NULL;
 205   AliAODTrack * track2 = NULL;
 206
 207
 208   for(int i = 0; i < tracks->GetEntriesFast(); i++) {
 209
 210     AliAODTrack * track = (AliAODTrack*)tracks->At(i);
 211     if (track->GetID() == l1) {
 212       track1 = track;
 213     } else if (track->GetID() == l2) {
 214       track2 = track;
 215     }
 216
 217     if(track1 && track2) break;
 218   }
 219
 220
 221   if(!track1 || !track2) return tag;
 222   if(track1->GetLabel() < 0 || track2->GetLabel() < 0) {
 223     //cout << "error balla"<< endl;
 224
 225   } else {
 226
 227     AliAODMCParticle * mcPart1 = dynamic_cast<AliAODMCParticle*>(arrayMC->At(track1->GetLabel()));
 228     AliAODMCParticle * mcPart2 = dynamic_cast<AliAODMCParticle*>(arrayMC->At(track2->GetLabel()));
 229
 230     if (mcPart1 && mcPart2) {
 231
 232       if(mcPart1->GetMother() == mcPart2->GetMother()) {
 233
 234         AliAODMCParticle * photon = dynamic_cast<AliAODMCParticle*>(arrayMC->At(mcPart1->GetMother()));
 235         Int_t motherIndex = photon->GetMother();
 236         tag= fAnaUtils->CheckOriginInAOD(&motherIndex, 1, arrayMC);
 237
 238         AliAODMCParticle * mother = dynamic_cast<AliAODMCParticle*>(arrayMC->At(motherIndex));
 239
 240
 241         if (photon->GetPdgCode() == 22 ) {
 242
 243           if (TMath::Abs(mother->GetPdgCode()) < 22 && TMath::Abs(mother->GetPdgCode()) != 11) {
 244
 245             fAnaUtils->SetTagBit(tag, AliMCAnalysisUtils::kMCPhoton);
 246             fAnaUtils->SetTagBit(tag, AliMCAnalysisUtils::kMCPrompt);
 247
 248             cout  <<"KKKKKKKKK "<< mother->GetPdgCode() << " " << mother->GetStatus() << " mi " << motherIndex << " daught: " << mother->GetNDaughters() << endl;
 249             cout << photon->GetStatus() << " " << photon->Pt() << " headerpt:" << mcHeader->GetPtHard() << " " <<mcPart1->GetMother() << endl;
 250
 251           }
 252         } else {
 253           cout << "FAKE "<< photon->PdgCode() << " " << photon->Pt() << " " <<  photon->GetStatus() << endl;
 254         }
 255
 256
 257
 258
 259
 260
 261         //Int_t parentId = mother->GetMother();
 262         //AliAODMCParticle * gp = dynamic_cast<AliAODMCParticle*>(arrayMC->At(parentId));
 263
 264
 265
 266
 267         if(! mother->IsPrimary()) {
 268           //cout << mother->GetPdgCode() << " " << mother->GetStatus() <<endl;
 269           //      cout << "other one? " << mother->PdgCode() << " " << gp->PdgCode() << " " << gp->GetStatus() << endl;
 270           //    cout << "yeay3  " <<  mother->GetPdgCode() << endl;
 271
 272         } else { //if (mother->IsPhysicalPrimary() ){
 273
 274
 275           if(mother->GetStatus() < 10) {
 276             //cout << "yeay4  " <<  mother->GetPdgCode() << endl;
 277             // cout << "pion? " << gp->PdgCode() << " " << gp->GetStatus() << endl;
 278             //cout << "???" << mother->PdgCode() << " " << mother->GetStatus() << " " << gp->PdgCode() << " " << gp->GetStatus() << endl;
 279           } else {
 280             //cout << "yeay5  " <<  mother->GetPdgCode() << endl;
 281             //cout << "other? " << mother->PdgCode() << " " << mother->GetStatus() << " " << gp->PdgCode() << " " << gp->GetStatus() << endl;
 282           }
 283
 284
 285           // } else {
 286           //   cout << "MI: " << mother->GetPdgCode() << " " << mother->GetStatus() <<endl;
 287         }
 288
 289       }
 290     }
 291   }
 292   //cout << "REturn tag " << tag << endl;
 293   return tag;
 294 }
 295
 296
 297 ///__________________________________________________________________________________
 298 AliAODPWG4ParticleCorrelation * AliAnalysisTaskGCPartToPWG4Part::AddToAOD(AliGammaConversionAODObject * aodO, TClonesArray * branch, TString detector) {
 299   new((*branch)[branch->GetEntriesFast()]) AliAODPWG4ParticleCorrelation(aodO->Px(), aodO->Py(), aodO->Pz(), aodO->E());
 300   AliAODPWG4ParticleCorrelation * photon = dynamic_cast<AliAODPWG4ParticleCorrelation*>(branch->Last());
 301   photon->SetTagged(aodO->IsTagged());
 302   photon->SetTrackLabel(aodO->GetLabel1(), aodO->GetLabel2());
 303   photon->SetDetector(detector);
 304   return photon;
 305 }
 306
 307 ///__________________________________________________________________________________
 308 AliAODPWG4ParticleCorrelation * AliAnalysisTaskGCPartToPWG4Part::AddToAOD(AliAODConversionParticle * aodO, TClonesArray * branch, TString detector) {
 309   new((*branch)[branch->GetEntriesFast()]) AliAODPWG4ParticleCorrelation(aodO->Px(), aodO->Py(), aodO->Pz(), aodO->E());
 310   AliAODPWG4ParticleCorrelation * photon = dynamic_cast<AliAODPWG4ParticleCorrelation*>(branch->Last());
 311   photon->SetTrackLabel(aodO->GetLabel1(), aodO->GetLabel2());
 312   photon->SetDetector(detector);
 313   return photon;
 314 }
 315
 316
 317
 318
 319 //_____________________________________________________________________
 320 void AliAnalysisTaskGCPartToPWG4Part::Terminate(Option_t *) {
 321   // Draw result to the screen
 322   // Called once at the end of the query
 323 }
 324
 325 //_____________________________________________________________________
 326 AliAODEvent * AliAnalysisTaskGCPartToPWG4Part::GetAODEvent() {
 327   //Get the AOD event from whereever it might be
 328   AliAODEvent * aodEvent = dynamic_cast<AliAODEvent*>(InputEvent());
 329   if(!aodEvent) {
 330     aodEvent = AODEvent();
 331   }
 332
 333   return aodEvent;
 334
 335 }
 336
 337 //_____________________________________________________________________
 338 TClonesArray * AliAnalysisTaskGCPartToPWG4Part::GetConversionGammas(const AliAODEvent * aodEvent) {
 339
 340   //Get Conversion gamma branch of AOD. First try standard AOD
 341   TClonesArray * convGamma = dynamic_cast<TClonesArray*>(aodEvent->FindListObject(fAODBranchName.Data()));
 342
 343
 344   //If it's there, send it back
 345   if(convGamma)  return convGamma;
 346
 347
 348   //If AOD not in standard file have to locate it in delta AOD
 349   if( !(fDeltaAODFileName.Length() > 0)  ) return NULL;
 350
 351   AliAODHandler * aodHandler = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
 352   if(aodHandler) {
 353     AliAODExtension * gExt = dynamic_cast<AliAODExtension*>(aodHandler->GetExtensions()->FindObject(fDeltaAODFileName));
 354     if(gExt) {
 355       AliAODEvent * gcEvent = gExt->GetAOD();
 356       return dynamic_cast<TClonesArray*>(gcEvent->FindListObject("GammaConv_gamma"));
 357     }
 358   }
 359   return NULL;
 360 }
 361
 362 //_________________________________________________________________________
 363 void AliAnalysisTaskGCPartToPWG4Part::CleanUp() {
 364   fAODPWG4Particles->Delete();
 365 }
 366
 367