PWGGA/EMCALJetTasks/AliJetEmbeddingTask.cxx

   1 // $Id$
   2 //
   3 // Jet embedding task.
   4 //
   5 // Author: Salvatore Aiola, Constantin Loizides
   6
   7 #include <TClonesArray.h>
   8 #include <TLorentzVector.h>
   9 #include <TRandom3.h>
  10
  11 #include "AliAnalysisManager.h"
  12 #include "AliVEvent.h"
  13 #include "AliVCluster.h"
  14 #include "AliEMCALDigit.h"
  15 #include "AliEMCALRecPoint.h"
  16 #include "AliPicoTrack.h"
  17 #include "AliEMCALGeometry.h"
  18 #include "AliLog.h"
  19
  20 #include "AliJetEmbeddingTask.h"
  21
  22 ClassImp(AliJetEmbeddingTask)
  23
  24 //________________________________________________________________________
  25 AliJetEmbeddingTask::AliJetEmbeddingTask() :
  26   AliAnalysisTaskSE("AliJetEmbeddingTask"),
  27   fGeomName(),
  28   fTracksName(),
  29   fOutTracksName(),
  30   fCaloName(),
  31   fOutCaloName(),
  32   fEtaMin(-1),
  33   fEtaMax(1),
  34   fPhiMin(0),
  35   fPhiMax(TMath::Pi() * 2),
  36   fPtMin(0),
  37   fPtMax(0),
  38   fCopyArray(kTRUE),
  39   fNEmbClusters(0),
  40   fNEmbTracks(0),
  41   fGeom(0),
  42   fClusters(0),
  43   fOutClusters(0),
  44   fTracks(0),
  45   fOutTracks(0)
  46 {
  47   // Default constructor.
  48 }
  49
  50 //________________________________________________________________________
  51 AliJetEmbeddingTask::AliJetEmbeddingTask(const char *name) :
  52   AliAnalysisTaskSE(name),
  53   fGeomName(""),
  54   fTracksName("PicoTracks"),
  55   fOutTracksName("PicoTracksEmbedded"),
  56   fCaloName("CaloClustersCorr"),
  57   fOutCaloName("CaloClustersCorrEmbedded"),
  58   fEtaMin(-1),
  59   fEtaMax(1),
  60   fPhiMin(0),
  61   fPhiMax(TMath::Pi() * 2),
  62   fPtMin(50),
  63   fPtMax(60),
  64   fCopyArray(kTRUE),
  65   fNEmbClusters(0),
  66   fNEmbTracks(1),
  67   fGeom(0),
  68   fClusters(0),
  69   fOutClusters(0),
  70   fTracks(0),
  71   fOutTracks(0)
  72 {
  73   // Standard constructor.
  74
  75 }
  76
  77 //________________________________________________________________________
  78 AliJetEmbeddingTask::~AliJetEmbeddingTask()
  79 {
  80   // Destructor
  81 }
  82
  83 //________________________________________________________________________
  84 void AliJetEmbeddingTask::Init()
  85 {
  86   // Init task.
  87
  88   if (fNEmbTracks > 0) {
  89     fTracks = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTracksName));
  90
  91     if (strcmp(fTracks->GetClass()->GetName(), "AliPicoTrack")) {
  92       AliError("Can only embed PicoTracks!");
  93       return;
  94     }
  95
  96     if (!fTracks) {
  97       AliError(Form("Couldn't retrieve tracks with name %s!", fTracksName.Data()));
  98       return;
  99     }
 100
 101     if (!fOutTracks) {
 102       fOutTracksName = fTracksName;
 103       if (fCopyArray) {
 104         fOutTracksName += "Embedded";
 105         fOutTracks = new TClonesArray(*fTracks);
 106         fOutTracks->SetName(fOutTracksName);
 107       }
 108       else {
 109         fOutTracks = fTracks;
 110       }
 111     }
 112
 113     if (fCopyArray) {
 114       if (!(InputEvent()->FindListObject(fOutTracksName)))
 115         InputEvent()->AddObject(fOutTracks);
 116       fOutTracks->Clear();
 117     }
 118   }
 119
 120   if (fNEmbClusters > 0) {
 121     fClusters = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fCaloName));
 122
 123     if (!fClusters) {
 124       AliError(Form("Couldn't retrieve clusters with name %s!", fCaloName.Data()));
 125       return;
 126     }
 127
 128     if (!fOutClusters) {
 129       fOutCaloName = fCaloName;
 130       if (fCopyArray) {
 131         fOutCaloName += "Embedded";
 132         fOutClusters = new TClonesArray(*fClusters);
 133         fOutClusters->SetName(fOutCaloName);
 134       }
 135       else {
 136         fOutClusters = fClusters;
 137       }
 138     }
 139
 140     if (fCopyArray) {
 141       if (!(InputEvent()->FindListObject(fOutCaloName)))
 142         InputEvent()->AddObject(fOutClusters);
 143       fOutClusters->Clear();
 144     }
 145   }
 146
 147     if (!fGeom) {
 148       if (fGeomName.Length()>0) {
 149         fGeom = AliEMCALGeometry::GetInstance(fGeomName);
 150         if (!fGeom)
 151           AliError(Form("Could not get geometry with name %s!", fGeomName.Data()));
 152       } else {
 153         fGeom = AliEMCALGeometry::GetInstance();
 154         if (!fGeom)
 155           AliError("Could not get default geometry!");
 156       }
 157     }
 158 }
 159
 160 //________________________________________________________________________
 161 void AliJetEmbeddingTask::Embed()
 162 {
 163   // Embed particles.
 164
 165   if (fNEmbClusters > 0 && fOutClusters) {
 166     const Int_t nClusters = fOutClusters->GetEntriesFast();
 167     TClonesArray digits("AliEMCALDigit", 1);
 168     for (Int_t i = 0; i < fNEmbClusters; ++i) {
 169       Double_t pt  = gRandom->Rndm() * (fPtMax - fPtMin) + fPtMin;
 170       Double_t eta = gRandom->Rndm() * (fEtaMax - fEtaMin) + fEtaMin;
 171       Double_t phi = gRandom->Rndm() * (fPhiMax - fPhiMin) + fPhiMin;
 172
 173       TLorentzVector nPart;
 174       nPart.SetPtEtaPhiM(pt, eta, phi, 0);
 175       Double_t e = nPart.E();
 176
 177       Int_t absId = 0;
 178       fGeom->GetAbsCellIdFromEtaPhi(eta, phi, absId);
 179
 180       AliEMCALDigit *digit = static_cast<AliEMCALDigit*>(digits.New(0));
 181       digit->SetId(absId);
 182       digit->SetIndexInList(0);
 183       digit->SetType(AliEMCALDigit::kHG);
 184       digit->SetAmplitude(e);
 185
 186       AliEMCALRecPoint recPoint;
 187       recPoint.AddDigit(*digit, e, kFALSE);
 188       recPoint.EvalGlobalPosition(0, &digits);
 189
 190       TVector3 gpos;
 191       recPoint.GetGlobalPosition(gpos);
 192       Float_t g[3];
 193       gpos.GetXYZ(g);
 194
 195       AliVCluster *cluster = static_cast<AliVCluster*>(fOutClusters->New(nClusters + i));
 196       cluster->SetType(AliVCluster::kEMCALClusterv1);
 197       cluster->SetE(recPoint.GetEnergy());
 198       cluster->SetPosition(g);
 199       cluster->SetNCells(1);
 200       UShort_t shortAbsId = absId;
 201       cluster->SetCellsAbsId(&shortAbsId);
 202       cluster->SetCellsAmplitudeFraction(NULL);
 203       cluster->SetID(nClusters + i);
 204       cluster->SetEmcCpvDistance(-1);
 205       cluster->SetChi2(100); // MC flag!
 206     }
 207   }
 208
 209   if (fNEmbTracks > 0 && fOutTracks) {
 210     Int_t nTracks = fOutTracks->GetEntriesFast();
 211     for (Int_t i = 0; i < fNEmbTracks; ++i) {
 212       Double_t pt  = gRandom->Rndm() * (fPtMax - fPtMin) + fPtMin;
 213       Double_t eta = gRandom->Rndm() * (fEtaMax - fEtaMin) + fEtaMin;
 214       Double_t phi = gRandom->Rndm() * (fPhiMax - fPhiMin) + fPhiMin;
 215
 216       new ((*fOutTracks)[nTracks + i]) AliPicoTrack(pt,
 217                                                     eta,
 218                                                     phi,
 219                                                     1,
 220                                                     100,    // MC flag!
 221                                                     0,
 222                                                     0,
 223                                                     kFALSE);
 224     }
 225   }
 226 }
 227
 228 //________________________________________________________________________
 229 void AliJetEmbeddingTask::UserExec(Option_t *)
 230 {
 231   // Execute per event.
 232
 233   Init();
 234
 235   Embed();
 236 }