PWGJE/EMCALJetTasks/AliJetConstituentTagCopier.cxx

   1 // $Id: AliJetConstituentTagCopier.cxx  $
   2 //
   3 // Copy tags from particle level constituent to detector level
   4 //
   5 // Author: S. Aiola
   6
   7 #include "AliJetConstituentTagCopier.h"
   8
   9 #include <TClonesArray.h>
  10 #include <TMath.h>
  11 #include <TLorentzVector.h>
  12
  13 #include "AliNamedArrayI.h"
  14 #include "AliVCluster.h"
  15 #include "AliVParticle.h"
  16 #include "AliEmcalParticle.h"
  17 #include "AliLog.h"
  18
  19 ClassImp(AliJetConstituentTagCopier)
  20
  21 //________________________________________________________________________
  22 AliJetConstituentTagCopier::AliJetConstituentTagCopier() :
  23   AliAnalysisTaskEmcal("AliJetConstituentTagCopier", kFALSE),
  24   fMCParticlesName(),
  25   fMCParticles(0),
  26   fMCParticlesMap(0)
  27 {
  28   // Default constructor.
  29 }
  30
  31 //________________________________________________________________________
  32 AliJetConstituentTagCopier::AliJetConstituentTagCopier(const char *name) :
  33   AliAnalysisTaskEmcal(name, kFALSE),
  34   fMCParticlesName("MCParticles"),
  35   fMCParticles(0),
  36   fMCParticlesMap(0)
  37 {
  38   // Standard constructor.
  39 }
  40
  41 //________________________________________________________________________
  42 AliJetConstituentTagCopier::~AliJetConstituentTagCopier()
  43 {
  44   // Destructor
  45 }
  46
  47 //________________________________________________________________________
  48 void AliJetConstituentTagCopier::ExecOnce()
  49 {
  50   // Execute once.
  51
  52   if (!fMCParticles) {
  53     fMCParticles = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fMCParticlesName));
  54     if (!fMCParticles) {
  55       AliError(Form("%s: Could not retrieve MC particles %s!", GetName(), fMCParticlesName.Data()));
  56       return;
  57     }
  58     else if (!fMCParticles->GetClass()->GetBaseClass("AliVParticle")) {
  59       AliError(Form("%s: Collection %s does not contain AliVParticle objects!", GetName(), fMCParticlesName.Data()));
  60       fMCParticles = 0;
  61       return;
  62     }
  63   }
  64
  65   if (!fMCParticlesMap) {
  66     fMCParticlesMap = dynamic_cast<AliNamedArrayI*>(InputEvent()->FindListObject(fMCParticlesName + "_Map"));
  67     // this is needed to map the MC labels with the indexes of the MC particle collection
  68     // if teh map is not given, the MC labels are assumed to be consistent with the indexes (which is not the case if AliEmcalMCTrackSelector is used)
  69     if (!fMCParticlesMap) {
  70       AliWarning(Form("%s: Could not retrieve map for MC particles %s! Will assume MC labels consistent with indexes...", GetName(), fMCParticlesName.Data()));
  71       fMCParticlesMap = new AliNamedArrayI("tracksMap",9999);
  72       for (Int_t i = 0; i < 9999; i++) {
  73         fMCParticlesMap->AddAt(i,i);
  74       }
  75     }
  76   }
  77
  78   AliAnalysisTaskEmcal::ExecOnce();
  79 }
  80
  81 //________________________________________________________________________
  82 Bool_t AliJetConstituentTagCopier::Run()
  83 {
  84   if (fTracks) {
  85     if (fTracks->GetClass()->GetBaseClass("AliVParticle"))
  86       DoTrackLoop(fTracks);
  87     else if (fTracks->GetClass()->GetBaseClass("AliEmcalParticle"))
  88       DoEmcalParticleLoop(fTracks);
  89     else
  90       AliError(Form("%s: Object type not recognized in collection %s. Nothing will be done.", GetName(), fTracks->GetName()));
  91   }
  92
  93   if (fCaloClusters) {
  94     if (fCaloClusters->GetClass()->GetBaseClass("AliVCluster"))
  95       DoClusterLoop(fCaloClusters);
  96     else if (fCaloClusters->GetClass()->GetBaseClass("AliEmcalParticle"))
  97       DoEmcalParticleLoop(fCaloClusters);
  98     else
  99       AliError(Form("%s: Object type not recognized in collection %s. Nothing will be done.", GetName(), fCaloClusters->GetName()));
 100   }
 101
 102   return kTRUE;
 103 }
 104
 105 //________________________________________________________________________
 106 void AliJetConstituentTagCopier::DoClusterLoop(TClonesArray *array)
 107 {
 108   Double_t totalEnergy = 0;
 109   for (Int_t i = 0; i < array->GetEntries(); i++) {
 110     AliVCluster *cluster = static_cast<AliVCluster*>(array->At(i));
 111     if (!cluster) {
 112       AliError(Form("%s: Could not get cluster %d", GetName(), i));
 113       continue;
 114     }
 115     if (!AcceptCluster(cluster))
 116       continue;
 117     Int_t mcLabel = cluster->GetLabel();
 118     if (mcLabel > 0) {
 119       TLorentzVector vect;
 120       cluster->GetMomentum(vect, fVertex);
 121       AliDebug(2, Form("Cluster %d, pt = %f, eta = %f, phi = %f, label = %d",
 122                        i, cluster->E(), vect.Eta(), vect.Phi(), mcLabel));
 123       totalEnergy += cluster->E();
 124       Int_t index = -1;
 125       if (mcLabel < fMCParticlesMap->GetSize())
 126         index = fMCParticlesMap->At(mcLabel);
 127       if (index < 0)
 128         continue;
 129       AliVParticle *part = static_cast<AliVParticle*>(fMCParticles->At(index));
 130       if (!part) {
 131         AliError(Form("%s: Could not get MC particle %d", GetName(), index));
 132         continue;
 133       }
 134       AliDebug(2, Form("Matched with particle %d, pt = %f, eta = %f, phi = %f",
 135                        index, part->E(), part->Eta(), part->Phi()));
 136       UInt_t bits = (UInt_t)part->TestBits(TObject::kBitMask);
 137       cluster->SetBit(bits);
 138     }
 139   }
 140
 141   AliDebug(2, Form("Total energy of MC clusters = %f",
 142                    totalEnergy));
 143 }
 144
 145 //________________________________________________________________________
 146 void AliJetConstituentTagCopier::DoTrackLoop(TClonesArray *array)
 147 {
 148   for (Int_t i = 0; i < array->GetEntries(); i++) {
 149     AliVParticle *track = static_cast<AliVParticle*>(array->At(i));
 150     if (!track) {
 151       AliError(Form("%s: Could not get track %d", GetName(), i));
 152       continue;
 153     }
 154     if (!AcceptTrack(track))
 155       continue;
 156     Int_t mcLabel = TMath::Abs(track->GetLabel());
 157     if (mcLabel != 0) {
 158       Int_t index = -1;
 159       if (mcLabel < fMCParticlesMap->GetSize())
 160         index = fMCParticlesMap->At(mcLabel);
 161       if (index < 0)
 162         continue;
 163       AliVParticle *part = static_cast<AliVParticle*>(fMCParticles->At(index));
 164       if (!part) {
 165         AliError(Form("%s: Could not get MC particle %d", GetName(), index));
 166         continue;
 167       }
 168       AliDebug(3, Form("Track %d, pt = %f, eta = %f, phi = %f, label = %d is matched with particle %d, pt = %f, eta = %f, phi = %f",
 169                        i, track->Pt(), track->Eta(), track->Phi(), mcLabel, index, part->Pt(), part->Eta(), part->Phi()));
 170       UInt_t bits = (UInt_t)part->TestBits(TObject::kBitMask);
 171       track->SetBit(bits);
 172     }
 173   }
 174 }
 175
 176 //________________________________________________________________________
 177 void AliJetConstituentTagCopier::DoEmcalParticleLoop(TClonesArray *array)
 178 {
 179   for (Int_t i = 0; i < array->GetEntries(); i++) {
 180     AliEmcalParticle *emcpart = static_cast<AliEmcalParticle*>(array->At(i));
 181     if (!emcpart) {
 182       AliError(Form("%s: Could not get EmcalParticle %d", GetName(), i));
 183       continue;
 184     }
 185     if (!AcceptEmcalPart(emcpart))
 186       continue;
 187     AliVCluster *cluster = emcpart->GetCluster();
 188     AliVParticle *track = emcpart->GetTrack();
 189     Int_t mcLabel = 0;
 190     if (cluster)
 191       mcLabel = cluster->GetLabel();
 192     else if (track)
 193       mcLabel = TMath::Abs(track->GetLabel());
 194     if (mcLabel != 0) {
 195       Int_t index = -1;
 196       if (mcLabel < fMCParticlesMap->GetSize())
 197         index = fMCParticlesMap->At(mcLabel);
 198       if (index < 0)
 199         continue;
 200       AliVParticle *part = static_cast<AliVParticle*>(fMCParticles->At(index));
 201       if (!part) {
 202         AliError(Form("%s: Could not get MC particle %d", GetName(), index));
 203         continue;
 204       }
 205       UInt_t bits = (UInt_t)part->TestBits(TObject::kBitMask);
 206       emcpart->SetBit(bits);
 207     }
 208   }
 209 }