AliPhysics  v5-06-21-01 (1eac791)
 All Classes Namespaces Files Functions Variables Enumerations Enumerator Macros
AliReducedHighPtEventCreator.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2015, 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 #include <map>
16 #include <vector>
17 
18 #include <TArrayD.h>
19 #include <TArrayI.h>
20 #include <TClonesArray.h>
21 #include <TMath.h>
22 #include <TObjArray.h>
23 #include <TString.h>
24 #include <TTree.h>
25 
26 #include "AliAnalysisUtils.h"
27 #include "AliAODEvent.h"
28 #include "AliAODMCParticle.h"
29 #include "AliAODTrack.h"
30 #include "AliCentrality.h"
31 #include "AliESDEvent.h"
32 #include "AliESDtrack.h"
33 #include "AliInputEventHandler.h"
34 #include "AliLog.h"
35 #include "AliMCEvent.h"
36 #include "AliPicoTrack.h"
37 #include "AliStack.h"
38 #include "AliVCluster.h"
39 #include "AliVTrack.h"
40 #include "AliVVertex.h"
41 #include "AliEmcalTriggerPatchInfo.h"
42 
44 #include "AliReducedEmcalCluster.h"
46 #include "AliReducedHighPtEvent.h"
47 #include "AliReducedMCHeader.h"
49 #include "AliReducedPatchInfo.h"
52 
54 ClassImp(HighPtTracks::AliReducedTrackSelectionContainer)
55 ClassImp(HighPtTracks::AliReducedHighPtEventCreator)
57 
58 namespace HighPtTracks {
59 
63 AliReducedHighPtEventCreator::AliReducedHighPtEventCreator():
64  AliAnalysisTaskEmcal(),
65  fOutputTree(NULL),
66  fOutputEvent(NULL),
67  fTrackSelections(NULL),
68  fSwapTriggerThresholds(kFALSE),
69  fMinClusterE(-1),
70  fMaxClusterE(1000),
71  fMinPt(-1),
72  fMaxPt(1000),
73  fMinEta(-1000),
74  fMaxEta(1000),
75  fCentralityMethod("V0A"),
76  fMinBiasSelection(AliVEvent::kINT7)
77 {
78 
79 }
80 
85 AliReducedHighPtEventCreator::AliReducedHighPtEventCreator(const char* name):
86  AliAnalysisTaskEmcal(name, kTRUE),
87  fOutputTree(NULL),
88  fOutputEvent(NULL),
89  fTrackSelections(NULL),
90  fSwapTriggerThresholds(kFALSE),
91  fMinClusterE(-1),
92  fMaxClusterE(1000),
93  fMinPt(-1),
94  fMaxPt(1000),
95  fMinEta(-1000),
96  fMaxEta(1000),
97  fCentralityMethod("V0A"),
98  fMinBiasSelection(AliVEvent::kINT7)
99 {
100  DefineOutput(2, TTree::Class());
101 
102  fTrackSelections = new TObjArray;
103  fTrackSelections->SetOwner(kTRUE);
104 }
105 
109 AliReducedHighPtEventCreator::~AliReducedHighPtEventCreator() {
110  if(fTrackSelections) delete fTrackSelections;
111 }
112 
116 void AliReducedHighPtEventCreator::UserCreateOutputObjects() {
117  AliAnalysisTaskEmcal::UserCreateOutputObjects();
118  OpenFile(2);
119 
120  fOutputTree = new TTree("ReducedEvent", "A simple reduced event");
121  fOutputEvent = new AliReducedHighPtEvent();
122  fOutputTree->Branch("ReducedEvent", "AliReducedHighPtEvent", fOutputEvent, 128000, 0);
123 
124  PostData(1, fOutput);
125  PostData(2, fOutputTree);
126 }
127 
138 Bool_t AliReducedHighPtEventCreator::Run() {
139  if(!fCaloClusters){
140  AliError("Cluster container missing");
141  return kFALSE;
142  }
143  if(!fTracks){
144  AliError("Track container missing");
145  return kFALSE;
146  }
147  if(!SelectEvent(fInputEvent)) return kFALSE;
148  new(fOutputEvent) AliReducedHighPtEvent(kTRUE);
149 
150  // Write event specific information
151  fOutputEvent->SetVertexZ(fInputEvent->GetPrimaryVertex()->GetZ());
152  AliCentrality *centralityHandler = fInputEvent->GetCentrality();
153  if(centralityHandler) fOutputEvent->SetCentralityPercentile(centralityHandler->GetCentralityPercentile(fCentralityMethod.Data()));
154  ConvertTriggerPatches(fTriggerPatchInfo, fOutputEvent->GetPatchContainer());
155  TString triggerString(fInputEvent->GetFiredTriggerClasses());
156  fOutputEvent->SetDecisionFromTriggerString(triggerString.Contains("EG1"), triggerString.Contains("EG2"), triggerString.Contains("EJ1"), triggerString.Contains("EJ2"));
157  fOutputEvent->SetMinBiasEvent(fInputHandler->IsEventSelected() & fMinBiasSelection);
158  fOutputEvent->SetRunNumber(fInputEvent->GetRunNumber());
159 
160  std::map<int, int> mcindexmap;
161  if(MCEvent()){
162  // Generate Monte-Carlo header
163  AliReducedMCHeader *mcheader = new AliReducedMCHeader();
164  if(fIsPythia){
165  mcheader->SetCrossSection(fXsection);
166  mcheader->SetNumberOfTrials(fNTrials);
167  mcheader->SetPtHard(fPtHard);
168  }
169  fOutputEvent->SetMonteCarloHeader(mcheader);
170 
171  // Convert Monte Carlo particles, fill index map
172  AliMCEvent *mcev = MCEvent();
173  int npart(0);
174  for(Int_t ipart = 0; ipart < mcev->GetNumberOfTracks(); ipart++){
175  AliVParticle *part = mcev->GetTrack(ipart);
176  Double_t pt(TMath::Abs(part->Pt())), eta(part->Eta());
177  if(pt < fMinPt || pt > fMaxPt) continue;
178  if(eta < fMinEta || eta > fMaxEta) continue;
179  if(!part->Charge()) continue;
180  if(part->IsA() == AliAODMCParticle::Class()){
181  AliAODMCParticle *aodpart = static_cast<AliAODMCParticle *>(part);
182  if(!aodpart->IsPhysicalPrimary()) continue;
183  } else {
184  // ESD - get information from stack
185  if(!mcev->Stack()->IsPhysicalPrimary(ipart)) continue;
186  }
187 
188  AliReducedGeneratedParticle *reducedgen = new AliReducedGeneratedParticle(npart, part->PdgCode(), part->Px(), part->Py(), part->Pz(), part->E());
189  fOutputEvent->AddReducedGeneratedParticle(reducedgen);
190  mcindexmap.insert(std::pair<int,int>(ipart, npart));
191  npart++;
192  }
193  }
194 
195  // Filter Clusters
196  std::map<int, int> clusterindexmap;
197  int ncluster(0);
198  Double_t vtxpos[3];fInputEvent->GetPrimaryVertex()->GetXYZ(vtxpos);
199  for(TIter cliter = TIter(fCaloClusters).Begin(); cliter != TIter::End(); ++cliter){
200  AliVCluster *incluster = static_cast<AliVCluster *>(*cliter);
201  if(!SelectCluster(incluster)) continue;
202  TLorentzVector clustervec;
203  incluster->GetMomentum(clustervec, vtxpos);
204  AliReducedEmcalCluster *redcluster = new AliReducedEmcalCluster(ncluster, incluster->E(), clustervec.Eta(), clustervec.Phi(), incluster->GetM02(), incluster->GetM20());
205  // Get leading 3 cell energies
206  TArrayD cellEnergies;
207  GetCellEnergies(incluster, cellEnergies);
208  TArrayI indices(cellEnergies.GetSize());
209  TMath::Sort(cellEnergies.GetSize(), cellEnergies.GetArray(), indices.GetArray(), kTRUE);
210  redcluster->SetLeadingCellEnergies(
211  cellEnergies[indices[0]],
212  cellEnergies.GetSize() > 1 ? cellEnergies[indices[1]] : 0,
213  cellEnergies.GetSize() > 2 ? cellEnergies[indices[2]] : 0
214  );
215  // Assing MC particles
216  if(MCEvent()){
217  for(Int_t ilab = 0; ilab < incluster->GetNLabels(); ilab++){
218  AliVParticle *assigned = MCEvent()->GetTrack(TMath::Abs(incluster->GetLabels()[ilab]));
219  if(!assigned) continue;
220  redcluster->AddTrueContributor(assigned->PdgCode(), assigned->Px(), assigned->Py(), assigned->Pz(), assigned->E());
221  }
222  }
223  fOutputEvent->AddReducedCluster(redcluster);
224  clusterindexmap.insert(std::pair<int,int>(incluster->GetID(), ncluster));
225  ncluster++;
226  }
227 
228  // Filter tracks
229  for(TIter trkiter = TIter(fTracks).Begin(); trkiter != TIter::End(); ++trkiter){
230  AliVTrack *trackRaw = static_cast<AliVTrack *>(*trkiter), *trackToCheck(NULL);
231  // handlinf for pico tracks as input
232  if(trackRaw->IsA() == AliPicoTrack::Class()){
233  AliPicoTrack *picotrack = static_cast<AliPicoTrack *>(trackRaw);
234  trackToCheck = picotrack->GetTrack();
235  } else {
236  trackToCheck = trackRaw;
237  }
238  FixTrackInputEvent(trackToCheck);
239  TArrayI cutIndices;
240  if(!SelectTrack(trackToCheck,cutIndices)) continue;
241  AliReducedReconstructedTrack *rectrack = new AliReducedReconstructedTrack;
242  double pvec[3];
243  trackToCheck->GetPxPyPz(pvec);
244  rectrack->SetMomentumVector(pvec[0], pvec[1], pvec[2]);
245  rectrack->SetCharge(static_cast<Char_t>(trackToCheck->Charge()));
246  for(Int_t icut = 0; icut < cutIndices.GetSize(); icut++) rectrack->SetTrackCuts(cutIndices[icut]);
247  if(MCEvent()){
248  // handle label
249  Int_t label = TMath::Abs(trackToCheck->GetLabel());
250  if(label < fMCEvent->GetNumberOfTracks()){
251  std::map<int,int>::iterator found = mcindexmap.find(label);
252  if(found != mcindexmap.end()){
253  rectrack->SetMatchedParticleIndex(found->second);
254  }
255  if(label >= 0) rectrack->SetGoodTrackLabel(kTRUE);
256  rectrack->SetTPCClusters(trackToCheck->GetTPCNcls());
257  rectrack->SetTPCCrossedRows(GetTPCCrossedRows(trackToCheck));
258  rectrack->SetTPCSharedClusters(trackToCheck->GetTPCSharedMapPtr()->CountBits());
259  rectrack->SetTPCFindableClusters(trackToCheck->GetTPCNclsF());
260  }
261  }
262  // handle clusters
263  Int_t clusterIndex = trackToCheck->GetEMCALcluster();
264  if(clusterIndex >= 0 && clusterIndex < fCaloClusters->GetEntries()){
265  std::map<int,int>::iterator found = clusterindexmap.find(clusterIndex);
266  if(found != clusterindexmap.end()){
267  rectrack->SetMatchedClusterIndex(found->second);
268  }
269  }
270  fOutputEvent->AddReducedReconstructedParticle(rectrack);
271  }
272 
273  fOutputTree->Fill();
274  PostData(2, fOutputTree);
275  return kTRUE;
276 }
277 
283 void AliReducedHighPtEventCreator::AddVirtualTrackSelection(
285  fTrackSelections->Add(new AliReducedTrackSelectionContainer(cutindex, sel));
286 }
287 
293 Bool_t AliReducedHighPtEventCreator::SelectEvent(AliVEvent* event) const {
294  const AliVVertex *primvtx = event->GetPrimaryVertex();
295  if(!primvtx || !primvtx->GetNContributors()) return kFALSE;
296  if(TMath::Abs(primvtx->GetZ()) > 10.) return kFALSE;
297  if(event->IsPileupFromSPD(3, 0.8, 3., 2., 5.)) return kFALSE;
298  AliAnalysisUtils eventSelUtil;
299  if(!eventSelUtil.IsVertexSelected2013pA(event)) return kFALSE;
300  return kTRUE;
301 }
302 
308 Bool_t AliReducedHighPtEventCreator::SelectCluster(const AliVCluster* clust) const {
309  if(!clust->IsEMCAL()) return kFALSE;
310  if(clust->E() < fMinClusterE || clust->E() > fMaxClusterE) return kFALSE;
311  return kTRUE;
312 }
313 
320 Int_t AliReducedHighPtEventCreator::SelectTrack(AliVTrack* track, TArrayI& cutindices) const {
321  std::vector<int> selected;
322 
323  double pt = TMath::Abs(track->Pt()), eta = track->Eta();
324  if(pt < fMinPt || pt > fMaxPt) return 0;
325  if(eta < fMinEta || eta > fMaxEta) return 0;
326 
327  for(TIter cutiter = TIter(fTrackSelections).Begin(); cutiter != TIter::End(); ++cutiter){
328  AliReducedTrackSelectionContainer *mycut = static_cast<AliReducedTrackSelectionContainer *>(*cutiter);
329  if(!mycut->GetTrackSelection()->IsTrackAccepted(track)) continue;
330  selected.push_back(mycut->GetIndex());
331  }
332  if(!selected.size()) return 0;
333  cutindices.Set(selected.size());
334  int counter = 0;
335  for(std::vector<int>::iterator inditer = selected.begin(); inditer != selected.end(); ++inditer){
336  cutindices[counter++] = *inditer;
337  }
338  return selected.size();
339 }
340 
346 Int_t AliReducedHighPtEventCreator::GetTPCCrossedRows(const AliVTrack* trk) const {
347  if(trk->IsA() == AliESDtrack::Class()){
348  return (static_cast<const AliESDtrack *>(trk))->GetTPCCrossedRows();
349  } else if(trk->IsA() == AliAODTrack::Class()){
350  return (static_cast<const AliAODTrack *>(trk))->GetTPCNCrossedRows();
351  }
352  return 0;
353 }
354 
360 void AliReducedHighPtEventCreator::GetCellEnergies(AliVCluster* emccluster, TArrayD& energies) const {
361  if(!fInputEvent->GetEMCALCells()) {
362  AliError("No EMCAL cells");
363  return;
364  }
365  AliDebug(2, Form("Number of cells: %d, array: %p", emccluster->GetNCells(), emccluster->GetCellsAbsId()));
366  energies.Set(emccluster->GetNCells());
367  for(int icell = 0; icell < emccluster->GetNCells(); icell++){
368  // printf("Cell ID: %d\n", emccluster->GetCellsAbsId()[icell]);
369  energies[icell] = fInputEvent->GetEMCALCells()->GetCellAmplitude(emccluster->GetCellsAbsId()[icell]);
370  }
371 }
372 
378 void AliReducedHighPtEventCreator::ConvertTriggerPatches(TClonesArray* patches,
379  AliReducedPatchContainer* cont) {
380  if(!patches){
381  AliError("Trigger patch container not found\n");
382  return;
383  }
384  for(TIter patchIter = TIter(patches).Begin(); patchIter != TIter::End(); ++patchIter){
385  AliEmcalTriggerPatchInfo *mypatch = static_cast<AliEmcalTriggerPatchInfo *>(*patchIter);
386  if(!mypatch->IsOfflineSimple() && mypatch->IsLevel0()) continue;
388  Bool_t isDefined(kFALSE);
389  if(mypatch->IsOfflineSimple()){
390  if(mypatch->IsGammaHighSimple()) { triggertype = AliReducedPatchContainer::kEMCGammaHigh; isDefined = kTRUE; }
391  if(mypatch->IsGammaLowSimple()) { triggertype = AliReducedPatchContainer::kEMCGammaLow; isDefined = kTRUE; }
392  if(mypatch->IsJetHighSimple()) { triggertype = AliReducedPatchContainer::kEMCJetHigh; isDefined = kTRUE; }
393  if(mypatch->IsJetLowSimple()) { triggertype = AliReducedPatchContainer::kEMCJetLow; isDefined = kTRUE; }
394  if(!isDefined){
395  AliDebug(2, "Unknown offline patch type");
396  continue;
397  }
398  AliDebug(2, Form("Adding offline patch of type %d", int(triggertype)));
399  cont->AddTriggerPatch(kTRUE, triggertype, mypatch->GetPatchE(), mypatch->GetADCAmp(), mypatch->GetEtaGeo(), mypatch->GetPhiGeo());
400  } else {
401  if(mypatch->IsGammaHigh()) { triggertype = fSwapTriggerThresholds ? AliReducedPatchContainer::kEMCGammaLow : AliReducedPatchContainer::kEMCGammaHigh; isDefined=kTRUE; }
402  if(mypatch->IsGammaLow()) { triggertype = fSwapTriggerThresholds ? AliReducedPatchContainer::kEMCGammaHigh : AliReducedPatchContainer::kEMCGammaLow; isDefined=kTRUE; }
403  if(mypatch->IsJetHigh()) { triggertype = fSwapTriggerThresholds ? AliReducedPatchContainer::kEMCJetLow : AliReducedPatchContainer::kEMCJetHigh; isDefined=kTRUE; }
404  if(mypatch->IsJetLow()) { triggertype = fSwapTriggerThresholds ? AliReducedPatchContainer::kEMCJetHigh : AliReducedPatchContainer::kEMCJetLow; isDefined=kTRUE; }
405  if(!isDefined){
406  AliDebug(2, "Unknown online patch type");
407  continue;
408  }
409  AliDebug(2, Form("Adding online patch of type %d", int(triggertype)));
410  cont->AddTriggerPatch(kFALSE, triggertype, mypatch->GetPatchE(), mypatch->GetADCAmp(), mypatch->GetEtaGeo(), mypatch->GetPhiGeo());
411  }
412  }
413 }
414 
420 void AliReducedHighPtEventCreator::FixTrackInputEvent(AliVTrack* trk) {
421  if(!trk->GetEvent()){
422  if(trk->IsA() == AliESDtrack::Class())
423  (static_cast<AliESDtrack *>(trk))->SetESDEvent(static_cast<AliESDEvent *>(fInputEvent));
424  else if(trk->IsA() == AliAODTrack::Class())
425  (static_cast<AliAODTrack *>(trk))->SetAODEvent(static_cast<AliAODEvent *>(fInputEvent));
426  else if(trk->IsA() == AliPicoTrack::Class()){
427  AliPicoTrack *mytrk = static_cast<AliPicoTrack *>(trk);
428  if(!mytrk->GetEvent()){
429  if(mytrk->GetTrack()->IsA() == AliESDtrack::Class())
430  (static_cast<AliESDtrack *>(mytrk->GetTrack()))->SetESDEvent(static_cast<AliESDEvent *>(fInputEvent));
431  else
432  (static_cast<AliAODTrack *>(mytrk->GetTrack()))->SetAODEvent(static_cast<AliAODEvent *>(fInputEvent));
433  }
434  }
435  }
436 }
437 
441 AliReducedTrackSelectionContainer::AliReducedTrackSelectionContainer():
442  fIndex(-1),
443  fTrackSelection(NULL)
444 {
445 }
446 
452 AliReducedTrackSelectionContainer::AliReducedTrackSelectionContainer(
454  fIndex(index),
455  fTrackSelection(sel)
456 {
457 }
458 
462 AliReducedTrackSelectionContainer::~AliReducedTrackSelectionContainer() {
463  if(fTrackSelection) delete fTrackSelection;
464 }
465 
466 } /* namespace HighPtTracks */
467 
Reduced event structure for high- analysis.
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
Declaration of class AliReducedHighPtEventCreator and AliReducedTrackSelectionContainer.
Declaration of class AlliReducedPatchInfo.
Declaration of class AliReducedPatchContainer, a container for reduced trigger patches.
Declaration of a reduced MC event header.
Declaration of class AliReducedGeneratedParticle.
Reduced information about reconstructed EMCAL clusters.
Declartion of class AliEMCalPtTaskVTrackSelection.